| Title: | Permutation Tests for Linear Models |
|---|---|
| Description: | Linear model functions using permutation tests. |
| Authors: | Bob Wheeler <[email protected]>, Marco Torchiano <[email protected]> |
| Maintainer: | Marco Torchiano <[email protected]> |
| License: | GPL (>= 2) |
| Version: | 2.1.1 |
| Built: | 2025-02-28 05:47:52 UTC |
| Source: | https://github.com/mtorchiano/lmperm |
Alfalfa data set
data(Alfalfa)data(Alfalfa)
Bob Wheeler [email protected]
aovp is aov modified to use permutation tests instead of normal
theory tests. Like aov, the ANOVA model is fitted by a call to lmp for each
stratum. Timing differences between aovp and aov are negligible.
aovp(formula, data = NULL, perm="Exact", seqs=FALSE, center=TRUE, projections = FALSE, qr = TRUE, contrasts = NULL, ...)aovp(formula, data = NULL, perm="Exact", seqs=FALSE, center=TRUE, projections = FALSE, qr = TRUE, contrasts = NULL, ...)
The arguments are the same as for aov.
Additional parameters may be included. They are described in the
"Additional Parameters" section below. These additional parameters
are the same as for lmp.
formula |
A formula specifying the model. |
data |
A data frame in which the variables specified in the formula will be found. If missing, the variables are searched for in the standard way. |
perm |
"Exact", "Prob", "SPR" will produce permutation probabilities. Anything else, such as "", will produce F-test probabilities. |
seqs |
If TRUE, will calculate sequential SS. If FALSE, unique SS. |
center |
If TRUE, numerical variables will be centered. |
projections |
Logical flag: should the projections be returned? |
qr |
Logical flag: should the QR decomposition be returned? |
contrasts |
A list of contrasts to be used for some of the factors
in the formula. These are not used for any |
... |
Arguments to be passed to |
The model Y=Xb+Zg+e is assumed, where X is the incidence matrix for fixed effects,
and Z is an incidence matrix for random effects, with columns representing the several error strata.
The aovp() algorithm projects Y into strata such that each stratum has a single error term, such
as a stratum defined by whole blocks. X is also projected so that the model
in this stratum becomes P(Y)=P(X)bi+ei.
The vector bi is divided into sources with dfj degrees of freedom for the
jth source, and summary(aovp()) will produce an ANOVA table for these sources
in the ith strata. See Venables and Ripley for details.
Either permutation test p-values or the usual F-test
p-values will be output. Polynomial model terms are collected into
sources, so that Y~A+B+I(A^2) will contain two sources, one for A with 2 df,
and one for B with 1 df. Sources for factors are treated as usual, and polynomial
terms and factors may be mixed in one model. The function poly.formula may
be used to create polynomial models, and the function multResp may be used to
create multiresponse matrices for the lhs from variables defined in data.
The Exact method will permute the values exactly. The Prob and SPR methods will approximate
the permutation distribution by randomly exchanging pairs of Y elements. The Exact
method will be used by default when the number of observations is less than
or equal to maxExact, otherwise Prob will be used.
Prob: Iterations terminate when the estimated standard error of the estimated
proportion p is less than p*Ca. The iteration continues until all sources and
coefficients meet this criterion or until maxIter is reached. See Anscome(1953)
for the origin of the criterion.
SPR: This method uses sequential probability ratio tests to decide between
the hypotheses p0 and p1 for a strength (alpha, beta) test. The test terminates
upon the acceptance or rejection of p0 or if maxIter is reached. See Wald (1947).
The power of the SPR is beta at p0 and increases to 1-beta at p1. Placing p0 and
p1 close together makes the cut off sharp.
Exact: This method generates all permutations of Y. It will generally be found
too time consuming for more than 10 or 11 observations, but note that aovp
may be used to divide the data into small enough blocks for which exact
permutation tests may be possible.
For Prob and SPR, one may set nCycle to unity to exchange all elements instead
of just pairs at each iteration, but there seems to be no advantage to doing this
unless the number of iterations is small – say less than 100.
The SS will be calculated sequentially, just as lm does; or they may be
calculated uniquely, which means that the SS for each source is calculated
conditionally on all other sources. This is SAS type III, which is also what drop1()
produces, except that drop1() will not drop main effects when interactions are present.
The parameter seqs may be used to override the default unique calculation behavior.
The usual output from aov, with permutation p-values instead of normal
theory p-values.
Iter |
For |
Accept |
For |
These are the same as for lmp.
If TRUE, settings such as sequential or unique will be printed. Default TRUE
If TRUE, SS/Resid SS will be used, otherwise SS. The default is TRUE
For Prob and SPR: The maximum number of iterations. Default 1000.
For Prob: Stop iterations when estimated standard error of the estimated p is less than Ca*p. Default 0.1
For SPR: Null hypothesis probability. Default 0.05
For SPR: Alternative hypothesis probability. Default 0.06
For SPR: Size of SPR test. Default 0.01
For SPR: Type II error for SPR test. Default 0.01
For Exact: maximum number of observations allowed. If data exceeds this, Prob is used. Default 10.
For Prob and SPR: Performs a complete random permutation, instead of pairwise exchanges, every nCycle cycles. Default 1000.
There is a vignette with more details and an example. To access it, type
vignette("lmPerm")
The default contrasts are set internally to (contr.sum, contr.poly), which means
that coefficients are either pairwise contrasts with the last level or polynomial contrasts.
Numerical variables should be centered in order to make them orthogonal to the constant when ANOVA is to be done.
Factors with many levels may cause problems for the methodology used by R. It is designed to work well with both balanced and unbalanced data, but it is best used for factors with no more than four or five levels. If difficulties occur, degrees of freedom for high order sources will be reduced, and sometimes the sources will be omitted entirely. An error message usually accompanies such behavior.
The variables inside the Error() term are treated in the order given. Thus
Error(A+B) will usually produce components in the order A,B, with B orthogonalized
with respect to A. This may cause confusion for unbalanced block structures.
This function will behave identically to aov() if the following parameters are set:
perm="", seq=TRUE, center=FALSE.
Bob Wheeler [email protected]
Statistical theory with engineering applications. Wiley. Table 7.4
Modern applied statistics with S-Plus, p176
Sequential analysis, Wiley, Sec. 5.3
## A simple randomized block example. # There are 7 blocks and 6 treatments. A first # analysis with blocks as a factor shows block to be significant and treatments not. data(Hald17.4) summary(aovp(Y~T+block,Hald17.4)) # Using the block to define a separate error strata tells a different story. summary(aovp(Y~T+Error(block),Hald17.4)) # There appears to be a linear trend in the blocks. This may be investigated by # extracting a linear component. The factor L was created by copying the block # factor and assigning it a linear contrast, like this # contrasts(L,1)<-contr.poly(7). The analysis then becomes. summary(aovp(Y~T+L+Error(block),Hald17.4)) # The L factor is not significant under permutation. It is significant when aov() # is used and the test is the usual F-test. ## From Venables and Ripley (2000) # This is a 2^3 factorial in the variables N,P,K. It is fractioned by using the # three way interaction, NPK, into two fractions of 4. Each of these fractions is # allocated to 3 blocks, making 6 blocks in all. An analysis with block as a # variable is the following. As may be seen, aovp() discards the confounded NPK interaction. data(NPK) summary(aovp(yield ~ block + N*P*K, NPK)) # Since the NPK interaction was confounded with blocks, the experimenter no doubt judged # it of lesser interest. It may however be examined by including blocks as an additional # error term as follows. The basic error level between blocks is of course larger than # that within blocks, so the NPK interaction would have to be substantially larger that # it would have had to be were it tested within blocks. summary(aovp(yield ~ N*P*K + Error(block), NPK)) # The SS calculated by aovp() are unique SS by default. That is, # they are sums of squares for the difference of a model with and without the source. The # resulting test is a test of the hypothesis that the source has no effect on the response. # Sequential SS, which are those produced by aov() may be obtained by setting the # parameter seqs=TRUE. simDesign is an unbalanced design created by the AlgDesign package. data(simDesign) summary(aovp(Y~.,simDesign)) summary(aovp(Y~.,simDesign,seqs=TRUE)) # Since there is only one stratum, these results are the same as would be obtained from anova(lmp(Y~.,simDesign)) # ANOVA for numerical variables. First using contrasts, then numeric variables. data(Federer276) summary(aovp(Plants~Variety*Treatment+Error(Rep/Plot),Federer276)) data(Federer276Numeric) summary(aovp(poly.formula(Plants~quad(Variety,Treatment)+Error(Rep/Plot)),Federer276Numeric)) # The coefficients and their p-values may be obtained by summaryC(aovp(Plants~Variety*Treatment+Error(Rep/Plot),Federer276))## A simple randomized block example. # There are 7 blocks and 6 treatments. A first # analysis with blocks as a factor shows block to be significant and treatments not. data(Hald17.4) summary(aovp(Y~T+block,Hald17.4)) # Using the block to define a separate error strata tells a different story. summary(aovp(Y~T+Error(block),Hald17.4)) # There appears to be a linear trend in the blocks. This may be investigated by # extracting a linear component. The factor L was created by copying the block # factor and assigning it a linear contrast, like this # contrasts(L,1)<-contr.poly(7). The analysis then becomes. summary(aovp(Y~T+L+Error(block),Hald17.4)) # The L factor is not significant under permutation. It is significant when aov() # is used and the test is the usual F-test. ## From Venables and Ripley (2000) # This is a 2^3 factorial in the variables N,P,K. It is fractioned by using the # three way interaction, NPK, into two fractions of 4. Each of these fractions is # allocated to 3 blocks, making 6 blocks in all. An analysis with block as a # variable is the following. As may be seen, aovp() discards the confounded NPK interaction. data(NPK) summary(aovp(yield ~ block + N*P*K, NPK)) # Since the NPK interaction was confounded with blocks, the experimenter no doubt judged # it of lesser interest. It may however be examined by including blocks as an additional # error term as follows. The basic error level between blocks is of course larger than # that within blocks, so the NPK interaction would have to be substantially larger that # it would have had to be were it tested within blocks. summary(aovp(yield ~ N*P*K + Error(block), NPK)) # The SS calculated by aovp() are unique SS by default. That is, # they are sums of squares for the difference of a model with and without the source. The # resulting test is a test of the hypothesis that the source has no effect on the response. # Sequential SS, which are those produced by aov() may be obtained by setting the # parameter seqs=TRUE. simDesign is an unbalanced design created by the AlgDesign package. data(simDesign) summary(aovp(Y~.,simDesign)) summary(aovp(Y~.,simDesign,seqs=TRUE)) # Since there is only one stratum, these results are the same as would be obtained from anova(lmp(Y~.,simDesign)) # ANOVA for numerical variables. First using contrasts, then numeric variables. data(Federer276) summary(aovp(Plants~Variety*Treatment+Error(Rep/Plot),Federer276)) data(Federer276Numeric) summary(aovp(poly.formula(Plants~quad(Variety,Treatment)+Error(Rep/Plot)),Federer276Numeric)) # The coefficients and their p-values may be obtained by summaryC(aovp(Plants~Variety*Treatment+Error(Rep/Plot),Federer276))
Factorial Experiment
data(CC164)data(CC164)
A three factor experiment involving lettuce plants. Response is the number of lettuce plants emerging after treatment by three levels of either nitrogen or phospate fertilizer.
Bob Wheeler [email protected]
Chochran, W. and Cox, G. (1957). Experimental Design, 2nd Ed. John Wiley & Sons, New York.
Data from a factorial experiment.
data(composite)data(composite)
A two variable factorial experiment to test the strength of a thermoplastic composite to laser intensity and tape speed.
Bob Wheeler [email protected]
Faraway, J. (2005). Linear Models with R. Chapman and Hall. N.Y.
Randomized block in three replicates
data(Federer276)data(Federer276)
A seed-gemination test on eight varieties of guayule. Four treatments were applied to lots of seed from the varieties.
A split plot design assigning varieties to whole plots was used. Treatments formed randomized blocks within whole plots. There were three replicates.
The whole plots were nested within replicates
Bob Wheeler [email protected]
Federer, W.T. (1955). Experimental Design, Macmillan, N.Y.
Randomized block in three replicates
data(Federer276)data(Federer276)
Varieties and Treatments are numeric variables.
A seed-gemination test on eight varieties of guayule. Four treatments were applied to lots of seed from the varieties.
A split plot design assigning varieties to whole plots was used. Treatments formed randomized blocks within whole plots. There were three replicates.
The whole plots were nested within replicates
Bob Wheeler [email protected]
Federer, W.T. (1955). Experimental Design, Macmillan, N.Y.
Gasoline blending mixture experiment
data(ghoctane)data(ghoctane)
An octane-blending experiment for a ternary system. Motor-octane ratings were obtained for a seven lattice point design, and two checkpoints.
Bob Wheeler [email protected]
Gorman, J.W. and Hinman, J.E. (1962). Simples lattice designs for multicomponent systems. Technometics. 4-4. 463-487. John Wiley & Sons, New York.
Data from Hald Table 17.4 on p 505.
data(Hald17.4)data(Hald17.4)
A randomized block experiment with 7 blocks and 6 treatments.
Bob Wheeler [email protected]
Hald, A. (1952). Statistical theory with engineering applications. Wiley. NY.
lmp is lm modified to use permutation tests instead of normal
theory tests. Like lm, it can be used to carry out regression,
single stratum analysis of variance and analysis of covariance . Timing
differences between lmp and lm are negligible.
lmp(formula, data, perm="Exact", seqs=FALSE, center=TRUE, subset, weights, na.action, method = "qr", model = TRUE, x = FALSE, y = FALSE, qr = TRUE, singular.ok = TRUE, contrasts = NULL, offset, ...)lmp(formula, data, perm="Exact", seqs=FALSE, center=TRUE, subset, weights, na.action, method = "qr", model = TRUE, x = FALSE, y = FALSE, qr = TRUE, singular.ok = TRUE, contrasts = NULL, offset, ...)
The arguments are mostly the same as for lm.
Additional parameters may be included. They are described in the
"Additional Parameters" section below. These additional parameters
are the same as for aovp.
formula |
a symbolic description of the model to be fit. |
data |
an optional data frame containing the variables
in the model. If not found in |
perm |
"Exact", "Prob", "SPR" will produce permutation probabilities. Anyting else, such as "", will produce F-test probabilites. |
seqs |
If TRUE, will calculate sequential SS. If FALSE, unique SS. |
center |
If TRUE will center numerical variables |
subset |
an optional vector specifying a subset of observations to be used in the fitting process. |
weights |
an optional vector of weights to be used
in the fitting process. If specified, weighted least squares is used
with weights |
na.action |
a function which indicates what should happen
when the data contain |
method |
the method to be used; for fitting, currently only
|
model, x, y, qr
|
logicals. If |
singular.ok |
logical. If |
contrasts |
an optional list. See the |
offset |
this can be used to specify an a priori
known component to be included in the linear predictor
during fitting. An |
... |
additional arguments to be passed. |
The usual regression model EY=Xb is assumed. The vector b is divided into sources
with dfi degrees of freedom for the ith source, and anova(lmp()) will produce an
ANOVA table for these sources. Either permutation test p-values or the usual F-test
p-values will be output. Polynomial model terms are collected into
sources, so that Y~A+B+I(A^2) will contain two sources, one for A with 2 df,
and one for B with 1 df. Sources for factors are treated as usual, and polynomial
terms and factors may be mixed in one model. The function poly.formula may
be used to create polynomial models, and the function multResp may
be used to create a multi-response matrix for the lhs from variables in data.
One may also use summary(lm()) to obtain coefficient estimates and
estimates of the permutation test p-values. The Exact method will permute
the values exactly. The Prob and SPR methods will approximate
the permutation distribution by randomly exchanging pairs of Y elements. The Exact
method will be used by default when the number of observations is less than
or equal to maxExact, otherwise Prob will be used.
Prob: Iterations terminate when the estimated standard error of the estimated
proportion p is less than p*Ca. The iteration continues until all sources and
coefficients meet this criterion or until maxIter is reached. See Anscome(1953)
for the origin of the criterion.
SPR: This method uses sequential probability ratio tests to decide between
the hypotheses p0 and p1 for a strength (alpha, beta) test. The test terminates
upon the acceptance or rejection of p0 or if maxIter is reached. See Wald (1947).
The power of the SPR is beta at p0 and increases to 1-beta at p1. Placing p0 and
p1 close together makes the cut off sharp.
Exact: This method generates all permutations of Y. It will generally be found
too time consuming for more than 10 or 11 observations, but note that aovp
may be used to divide the data into small enough blocks for which exact
permutation tests may be possible.
For Prob and SPR, one may set nCycle to unity to exchange all elements instead
of just pairs at each iteration, but there seems to be no advantage to doing this
unless the number of iterations is small – say less than 100.
The SS will be calculated sequentially, just as lm() does; or they may be
calculated uniquely, which means that the SS for each source is calculated
conditionally on all other sources. This is SAS type III, which is also what drop1()
produces, except that drop1() will not drop main effects when interactions are present.
The parameter seqs may be used to override the default unique calculation behavior.
The usual output from lm, with permutation p-values or F-test
p-values. The p-values for the coefficients are of necessity, two-sided.
Iter |
For Prob and SPR: The number of iterations until the criterion is met. |
Accept |
For SPR: 1 the p0 hypothesis is accepted, and 0 for rejection or when no decision before |
These are the same as for aovp.
If TRUE, settings such as sequential or unique will be printed. Default TRUE
If TRUE, SS/Resid SS will be used, otherwise SS. The default is TRUE
For Prob and SPR: The maximum number of iterations. Default 1000.
For Prob: Stop iterations when estimated standard error of the estimated p is less than Ca*p. Default 0.1
For SPR: Null hypothesis probability. Default 0.05
For SPR: Alternative hypothesis probability. Default 0.06
For SPR: Size of SPR test. Default 0.01
For SPR: Type II error for SPR test. Default 0.01
For Exact: maximum number of observations allowed. If data exceeds this, Prob is used. Default 10.
For Prob and SPR: Performs a complete random permutation, instead of pairwise exchanges, every nCycle cycles. Default 1000.
There is a vignette with more details and an example. To access it, type
vignette("lmPerm")
The default contrasts are set internally to (contr.sum, contr.poly), which means
that factor coefficients are either pairwise contrasts with the last level or polynomial contrasts.
Numerical variables should be centered in order to make them orthogonal to the constant when ANOVA is to be done.
This function will behave identically to lm() if the following parameters are set:
perm="", seq=TRUE, center=FALSE. An exception for multiple responses is that an
ANOVA table for each response is output instead of a call to anova.mlm().
Bob Wheeler [email protected]
Experimental Design, 2nd Ed.
John Wiley & Sons, New York.
Sequential analysis, Wiley, Sec. 5.3
"Product improvement by application of Taguchi methods." in American Supplier Institute News (special symposium ed.) Dearborn, MI. American Supplier Institute. 11-16.
Signal-to-noise ratios, performance criteria, and transformations. Technometics. 30-1. 1-17.
# 3x3 factorial with ordered factors, each is average of 12. # This is a saturated design with no df for error. The results tend to support # Cochran and Cox who used a guessed residual SS for their analysis. The design # is balanced, so the sequential SS are the same as the unique SS. data(CC164) summary(lmp(y ~ N * P, data = CC164, perm="")) # F-value output as if lm() was used. summary(lmp(y ~ N * P, data = CC164,)) # Default, using "Exact" if possible. summary(lmp(y ~ N * P, data = CC164, perm="SPR")) anova(lmp(y ~ N * P, data = CC164)) # A two level factorial. The artificial data is N(0,1) with an effect of # 1.5 added to factor X4. When the number of iterations are small, as in # this case, using nCycle=1 is advantageous. X<-expand.grid(X1=1:2,X2=1:2,X3=1:2,X4=1:2) X$Y<-c(0.99,1.34,0.88,1.94,0.63,0.29,-0.78,-0.89,0.43,-0.03,0.50,1.66,1.65,1.78,1.31,1.51) summary(lmp(Y~(X1+X2+X3+X4)^2,X,"SP")) # The prob method is used because "SP" is not recognized. summary(lmp(Y~(X1+X2+X3+X4)^2,X,"SPR")) summary(lmp(Y~(X1+X2+X3+X4)^2,X,"SPR",nCycle=1)) #An additional parameter being passed. # A saturated design with 15 variables in 16 runs. The orginal analysis by Quinlan pooled the mean # squares from the 7 smallest effcts and found many variables to be significant. Box, reanalyzed # the data using half-normal plots and found only variables E and G to be important. The permutation # analysis agrees with this conclusion. data(Quinlan) summary(lmp(SN~.,Quinlan)) # A design containing both a polynomial variable and a factor data(simDesignPartNumeric) anova(lmp(poly.formula(Y~quad(A,B)+C),simDesignPartNumeric))# 3x3 factorial with ordered factors, each is average of 12. # This is a saturated design with no df for error. The results tend to support # Cochran and Cox who used a guessed residual SS for their analysis. The design # is balanced, so the sequential SS are the same as the unique SS. data(CC164) summary(lmp(y ~ N * P, data = CC164, perm="")) # F-value output as if lm() was used. summary(lmp(y ~ N * P, data = CC164,)) # Default, using "Exact" if possible. summary(lmp(y ~ N * P, data = CC164, perm="SPR")) anova(lmp(y ~ N * P, data = CC164)) # A two level factorial. The artificial data is N(0,1) with an effect of # 1.5 added to factor X4. When the number of iterations are small, as in # this case, using nCycle=1 is advantageous. X<-expand.grid(X1=1:2,X2=1:2,X3=1:2,X4=1:2) X$Y<-c(0.99,1.34,0.88,1.94,0.63,0.29,-0.78,-0.89,0.43,-0.03,0.50,1.66,1.65,1.78,1.31,1.51) summary(lmp(Y~(X1+X2+X3+X4)^2,X,"SP")) # The prob method is used because "SP" is not recognized. summary(lmp(Y~(X1+X2+X3+X4)^2,X,"SPR")) summary(lmp(Y~(X1+X2+X3+X4)^2,X,"SPR",nCycle=1)) #An additional parameter being passed. # A saturated design with 15 variables in 16 runs. The orginal analysis by Quinlan pooled the mean # squares from the 7 smallest effcts and found many variables to be significant. Box, reanalyzed # the data using half-normal plots and found only variables E and G to be important. The permutation # analysis agrees with this conclusion. data(Quinlan) summary(lmp(SN~.,Quinlan)) # A design containing both a polynomial variable and a factor data(simDesignPartNumeric) anova(lmp(poly.formula(Y~quad(A,B)+C),simDesignPartNumeric))
Lizard data
data(manly126)data(manly126)
Two factor data showing ant consumption by lizards. In addition, three aditional columns, v1,v2,v3, are included representing the additions used by Manly to calculate power.
Bob Wheeler [email protected]
Manly, B.F.J. (1998). Randomization, bootstrap and Monte Carlo methods in Biology. 2nd Ed. Chapman & Hall, London.
Ephemeroptera counts from New Zealand streams.
data(manly136)data(manly136)
Bob Wheeler [email protected]
Manly, B.F. (1998). Randomization, bootstrap and monte carlo methods in biology. Chapman and Hall, NY.
This function creates a multiple response matrix for its argument variables. When used on the lhs of the formula in lmp() or aovp() it will create a matrix containing one or more response columns from variables defined in the data argument.
multResp(...)multResp(...)
... |
Variable names separated by commas. |
A matrix with named columns
Bob Wheeler [email protected]
Please cite this program as follows:
Wheeler, R.E. (2010). multResp() lmPerm. The R project for statistical computing http://www.r-project.org/
A<-1:5 B<-1:5 multResp(A,B) data(Plasma) anova(lmp(multResp(Amin,Pct,sinpoly)~.,Plasma))A<-1:5 B<-1:5 multResp(A,B) data(Plasma) anova(lmp(multResp(Amin,Pct,sinpoly)~.,Plasma))
Data from Venables and Ripley, Table 6.1, p176.
data(NPK)data(NPK)
A fractional factorial with blocks confounded with the three way interaction.
Bob Wheeler [email protected]
Venables, W.N. and Ripley, B.D. (2000) Modern applied statistics with S-Plus. 3rd ed. Springer, NY.
Generates permutations a pair at a time.
permute(N=4,K=1,initialize=0)permute(N=4,K=1,initialize=0)
N |
The number of elements in the permutation vector. The elements are numbered 1:N |
K |
The number of exchange pairs to be returned at each access |
initialize |
Set to 1 to initialize. Set to 0 for additional pairs. |
On first call, set initialize to 1. On subsequent accesses, initialize should be 0. The function returns at most K exchange pairs at each access. The function should be repeatedly called until result is false.
The function returns a list.
result |
Returns false when the last permutation has been generated or a failure has occurred. Returns true for success and at most K new pairs in vec. |
N |
The input parameter |
K |
The input parameter |
vec |
A 2*K vector of element numbers to be exchanged – permutations may be obtained by successively applying these pairs. |
initialize |
Always returns 1 |
count |
The number of pairs in vec. It may be less than K for the final set. |
This is the minimum time and effort routine. It is used for Exact permutation calculations. The present function is a wrapper for the C code and may be useful for other purposes. The code follows Reingole and Nievergelt (1977).
Bob Wheeler [email protected]
Combinatorial Algorithms Theory and Practice. Prentice Hall, New Jersey. p170.
An experiment designed to guide optimization of a nitride etch process
data(Plasma)data(Plasma)
An orthogonal arry, capable of supporting a linear model only. The factors are denoted by their units of measurement. They were W, the power applied to the cathode; mTorr, the pressure in the reaction chamber; cm, gap between anode and cathode; sscm, flow if the reactant gas. The three responses: Amin, etch rate; Pct, unifomity; sinpoly, selectivity.
Bob Wheeler [email protected]
Vardeman, Stephen B. (1994).Statistics for engineering problem solving. PWS publishing Co. Boston. p 596
Formulas are expanded to accommodate special functions for continuous and mixture variables.
poly.formula(frml)poly.formula(frml)
frml |
A formula using ~ in the usual way. |
This function expands formulas to accommodate polynomial models for which R has minimal support. Assuming for illustration that there are three variables, A, B, and C, the following expressions may be used.
All agruments to quad(), cubic(), and cubicS() must be numeric.
quad(A,B,C) makes
cubic(A,B,C) makes
cubicS(A,B,C) makes
The cubicS() function produces a non-singular representation of a cubic model, when the
variables are mixture variables, that is when the rows of data sum to a constant
value, usually 1.0. Because of the mixture constraint, models containing mixture variables
should not have a constant term. The linear and quadratic models for mixture variables
A, B, and C are given by and respectively. See Gorman and Hinman [1962] for
details.
An expanded formula is returned.
Bob Wheeler [email protected]
Please cite this program as follows:
Wheeler, R.E. (2010). poly.formula lmPerm. The R project for statistical computing http://www.r-project.org/
Gorman, J.W. and Hinman, J.E. (1962). Simplex lattice designs for multicomponent systems. Technometrics. 4-4. 463-487.
poly.formula(y~quad(A,B,C)+Error(block))poly.formula(y~quad(A,B,C)+Error(block))
Data from a Symposium on Taguchi Methods, analyzed by G.E.P. Box
data(Quinlan)data(Quinlan)
A saturated design: L16 in Taguchi termonology, but just a 2^(15-11)III design using the symbolism in Box (1978). Reanalyzed by Box (1988).
Bob Wheeler [email protected]
"Product improvement by application of Taguchi methods." in American Supplier Institute News (special symposium ed.) Dearborn, MI. American Supplier Institute. 11-16.
Statistics for Experimenters. Wiley, N.Y
Signal-to-noise ratios, performance criteria, and transformations. Technometics. 30-1. 1-17.
Data from Scheffe, p140
data(ratGenotype)data(ratGenotype)
An unbalanced design. The data originally appeared in Bailey (1953), and was reproduced in Scheffe (1959) with the description:
“...gives the weight of hybrid female rats in a foster-nursing experiment with four types of rats. (The weights are litter averges in grams at 28 days. The within litter variance was obviously negligible compared to the between-litter variance.) The factors in the two-way laout are the genotype of the foster mother and that of the litter. ”
Bob Wheeler [email protected]
The Analysis of Variance, Wiley, NY.
“The inheritance of Maternal Influences on the Growth of the Rat”. Ph.D. thesis, Univ. California
A three factor design created using AlgDesign, with a normal response.
data(simDesign)data(simDesign)
Bob Wheeler [email protected]
A three factor design created using AlgDesign, with a normal response.
data(simDesignPartNumeric)data(simDesignPartNumeric)
Bob Wheeler [email protected]
Replaces corresponding functions in base package.
## S3 method for class 'lmp' summary(object, correlation = FALSE, symbolic.cor = FALSE, ...) ## S3 method for class 'mlmp' summary(object, ...) ## S3 method for class 'summary.lmp' print(x, digits = max(3, getOption("digits") - 3), symbolic.cor = x$symbolic.cor, signif.stars= getOption("show.signif.stars"), ...) ## S3 method for class 'aovp' summary(object, intercept = FALSE, split, expand.split = TRUE, keep.zero.df = TRUE, ...) ## S3 method for class 'lmp' anova(object, ...)## S3 method for class 'lmp' summary(object, correlation = FALSE, symbolic.cor = FALSE, ...) ## S3 method for class 'mlmp' summary(object, ...) ## S3 method for class 'summary.lmp' print(x, digits = max(3, getOption("digits") - 3), symbolic.cor = x$symbolic.cor, signif.stars= getOption("show.signif.stars"), ...) ## S3 method for class 'aovp' summary(object, intercept = FALSE, split, expand.split = TRUE, keep.zero.df = TRUE, ...) ## S3 method for class 'lmp' anova(object, ...)
Same as for the corresponding functions in base package:
object |
an object of class |
x |
an object of class |
correlation |
logical; if |
digits |
the number of significant digits to use when printing. |
symbolic.cor |
logical. If |
signif.stars |
logical. If |
intercept |
logical: should intercept terms be included? |
split |
an optional named list, with names corresponding to terms in the model. Each component is itself a list with integer components giving contrasts whose contributions are to be summed. |
expand.split |
logical: should the split apply also to interactions involving the factor? |
keep.zero.df |
logical: should terms with no degrees of freedom be included? |
... |
further arguments passed to or from other methods. |
These modified functions are needed because the perm values, which are attached to the object, replace the usual test columns in the output from these functions.
Bob Wheeler [email protected]
Summarize an analysis of variance model showing coefficients
## S3 method for class 'aovlist': summaryC(object, ...)## S3 method for class 'aovlist': summaryC(object, ...)
object |
an object of class |
... |
further arguments passed on. |
Summarizes aovp output in terms of coefficients with p-values.
Bob Wheeler [email protected]
data(Alfalfa) summaryC(aovp(Yield~Variety*Date+Error(Block/Variety),Alfalfa))data(Alfalfa) summaryC(aovp(Yield~Variety*Date+Error(Block/Variety),Alfalfa))
A three factor experiment with non-normal response
data(wool)data(wool)
There are three, 3-level factors:
length of test specimen in mm
amplitude of loading cycle in mm
the load in g
The response was the number of cycles to failure of worsted yarn under cycles of repeated loading
Bob Wheeler [email protected]
Box, G.E.P. and Cox, D.R. (1964) An analysis of transformations (with discussion). J.R. Statists. Soc. B. 26. 211-246.