|
permute(npc)
|
R Documentation
|
Permute
Description
This function provides a Monte Carlo distribution of a linear test statistic.
Usage
permute(X,Y,A,B=1000)
Arguments
|
X
|
A vector/array of contrasts/weights of
dimension n x C (or
length = n if X is a vector), where n is the total number of
observations and C is the number of contrasts/weights.
|
|
Y
|
A vector/array of responses, with
dimension n x p (or
length = n if Y is a vector), where p is the number of variables.
|
|
A
|
A vector of
alternatives of length p. Alternative entries are: ‘1’
meaning ‘greater than’, ‘0’ meaning ‘not equal
to’ and ‘-1’
meaning ‘less than’.
|
|
B
|
Number of Monte
Carlo permutations to be done.
|
Details
Value
A matrix
T with the Monte Carlo distribution of the test
statistic whose dimensions are [B+1,
p, C] in the most general case
when p > 1 and C > 1; [B+1, p] when p > 1 and C
= 1; and [B+1, C] when p = 1 and C
> 1; if p = C = 1, T is a
vector of length B+1.
See Also
combine,t2p
## Example on producing plastic film from Krzanowski (1998, p. 381)
## Two sample problem with three (dependent) variables. Partial tests are required with one-sided alternatives.
n1<-9
n2<-11
n<-n1+n2
tear <- c(6.5, 6.2, 5.8, 6.5, 6.5, 6.9, 7.2, 6.9, 6.1, 6.3,
6.7, 6.6, 7.2, 7.1, 6.8, 7.1, 7.0, 7.2, 7.5, 7.6)
gloss <- c(9.5, 9.9, 9.6, 9.6, 9.2, 9.1, 10.0, 9.9, 9.5, 9.4,
9.1, 9.3, 8.3, 8.4, 8.5, 9.2, 8.8, 9.7, 10.1, 9.2)
opacity <- c(4.4, 6.4, 3.0, 4.1, 0.8, 5.7, 2.0, 3.9, 1.9, 5.7,
2.8, 4.1, 3.8, 1.6, 3.4, 8.4, 5.2, 6.9, 2.7, 1.9)
Y <- cbind(tear, gloss, opacity)
X <- c(rep(1/n1,n1),rep(-1/n2,n2))
# X: vector of contrasts: the test statistics is the difference of sampling means for all variables
A= c(1,1,-1)
# array specifying the alternatives (“greater” for ‘tear’ and ‘gloss’, “less” for ‘opacity’
source("permute.r")source("combine.r")source("t2p.r")
set.seed(3)
T<-permute(X,Y,A,B=1000)
dim(T)
# [1] 1001 3
P<-t2p(T)
P[1,]
# [1] 0.994 0.014 0.252 # partial p-values;
par(mfrow=c(1,p)) # data representation;
for(j in 1:p){boxplot(Y[,j]~X,main=colnames(Y)[j],names=c("Sample2","Sample1"))}
T1 <- combine(P,fun="Fisher",which=2,W=c(0.5,0.2,0.3))
dim(T1)
#[1] 1001 1
global.p <- t2p(T1)[1]
global.p
# [1] 0.248 # global p-value.
######### ONE WAY ANOVA:
######### WAIS score example, from Mack & Wolfe (1981).
n = 3 ; C = 5;
Y = c(8.62,9.94,10.06,9.85,10.43,11.31,9.98,10.69,11.4,9.12,9.89,10.57,4.8,9.18,9.27)
X = array(0,dim=c(n*C,C))
for(k in 1:C){X[(n*(k-1)+1):(n*k),k] = 1/n
}
# X: matrix whose columns are indicator vectors of each group; the partial tests are squared sample means
x11() ; boxplot(Y~rep(1:5,each=3))
set.seed(9)
T = permute(X,Y,A = 0,B=1000)
dim(T)
# [1] 1001 5
T[1,]
# [1] 91.0116 110.8809 114.2761 97.2196 60.0625 # squared group means
T1 = combine(T,fun="Direct",which=2,W=rep(n,C)) # T1 is equivalent to the between-group deviance
p.glob = t2p(T1)[1]
p.glob
# [1] 0.052 # global p-value
T1[1] - n*C*mean(Y)^2 # between-groups deviance
#[1] 16.55796
[Package npc
Index]