2014-02-14 44 views
0

我有一个data.frame与几个变量X1,X2 ...和一个分组变量“站点”我想找到X1比例站点== 1大于X1与现场== 2,我可以做到这一点与固定数目的位点的水平,并在每个时间变量,但我想泛指任何数目的水平和几个变量的,下面是一个例子:自动化数据框中的子集之间的两两比较R

# Generate data 
set.seed(20130226) 

n <- 100 
x1 <- matrix(c(rnorm(n, mean = 2),rnorm(n, mean = 5)),ncol=2) 
x2 <- matrix(c(rnorm(n, mean = 1), rnorm(n, mean = 4)),ncol=2) 
x3 <- matrix(c(rnorm(n, mean = 3), rnorm(n, mean = 3)),ncol=2) 
xx <- data.frame(x1,site=1) 
xx <- rbind(xx, data.frame(x2,site=2)) 
xx <- rbind(xx, data.frame(x3,site=3)) 

# comparisons 

s <- unique(xx$site) 
me1 <- with(xx,xx[site==s[1],]) 
me2<- with(xx,xx[site==s[2],]) 
me3<- with(xx,xx[site==s[3],]) 

Pg1.gt.g2 <- sum(me1[,c("X1")]>me2[,c("X1")])/nrow(me1) 
Pg1.gt.g3 <- sum(me1[,c("X1")]>me3[,c("X1")])/nrow(me1) 
Pg2.gt.g3 <- sum(me2[,c("X1")]>me3[,c("X1")])/nrow(me1) 

# build table 
comp1 <- data.frame(Group=c(paste(s[1],">",s[2]),paste(s[1],">",s[3]),paste(s[2],">",s[3])), P=c(Pg1.gt.g2, Pg1.gt.g3,Pg2.gt.g3)) 

print(comp1) 

我不知道如何做到这一点为不同数量的组和几个变量,也许使用plyr

谢谢!

回答

1

我的数据重塑成一个矩阵,其中每一列代表一组:

# Unique sites 
s <- unique(xx$site) 

# Columns are each group, data are X1 values 
mat <- do.call(cbind, lapply(split(xx, xx$site), function(x) x$X1)) 

# Compare all pairs of sites 
do.call(rbind, apply(combn(seq_along(s), 2), 2, 
        function(x) data.frame(g1=s[x[1]], g2=s[x[2]], 
              prop=sum(mat[,x[1]] > mat[,x[2]])/nrow(mat)))) 

# g1 g2 prop 
# 1 1 2 0.83 
# 2 1 3 0.20 
# 3 2 3 0.09 
+0

聪明,无法理解非常好最后do.call。无论如何,这是为X1,但如果你有一个未定义数量的变量? – Leosar

+0

如果你想为不同的变量做它,你需要改变你创建'mat'的方式 - 其他的都是一样的。如果你想计算许多不同变量的表格,你可以使这个代码成为一个函数,并传递变量名称(这里应该是“X1”)。 – josliber