如何使用测试统计分布对列名进行 1000 次排列？

Question

假设我有一个这样的矩阵

dat <- read.table(text = "   code.1 code.2 code.3 code.4
1     82     93     NA     NA
2     15     85     93     NA
3     93     89     NA     NA
4     81     NA     NA     NA",
                  header = TRUE, stringsAsFactors = FALSE)
dat2=data.matrix(dat)

实际上，我的矩阵有 132 列和大约 15000 行。 我的列名如下所示：NoD_14569_norm.1 NoD_14569_norm.2 NoD_14569_norm.3 NoD_14581_30mM.1 NoD_14581_30mM.2 NoD_14581_30mM.3

我想要做的是为我的列名创建 1000 个随机排列，其中矩阵中的所有内容都将保持不变，除非会有新的列名顺序。

例如，列名的一种排列/改组会给我这个：

  code.2 code.4 code.1 code.3
1     82     93     NA     NA
2     15     85     93     NA
3     93     89     NA     NA
4     81     NA     NA     NA

目标是对 1000 个数据帧中的每一个执行以下代码

subject="all_replicate"
targets<-readTargets(paste(PhenotypeDir,"hg_sg_",subject,"_target.txt", sep=''))
Treat <- factor(targets$Treatment,levels=c("C","T"))
Replicates <- factor(targets$rep)
design <- model.matrix(~Replicates+Treat)
corfit <- duplicateCorrelation(dat2, block = targets$Subject)
corfit$consensus.correlation
fit <-lmFit(dat2,design,block=targets$Subject,correlation=corfit$consensus.correlation)
fit<-eBayes(fit)
y1=topTable(fit, coef="TreatT", n=nrow(genes),adjust.method="BH",genelist=genes)

在 y1 内部有包含 p 值的列名 P.value，我想绘制所有上述 1000 个列名排列的分布。

请指教

Answer 1

列名的随机排序很简单：

set.seed(42)
# manyorders <- replicate(1000, sample(colnames(dat2)), simplify=FALSE)
# set.seed(42)
manyorders <- replicate(1000, sample(colnames(dat2)), simplify=FALSE)
head(manyorders)
# [[1]]
# [1] "code.4" "code.3" "code.1" "code.2"
# [[2]]
# [1] "code.3" "code.2" "code.4" "code.1"
# [[3]]
# [1] "code.3" "code.4" "code.1" "code.2"
# [[4]]
# [1] "code.4" "code.1" "code.3" "code.2"
# [[5]]
# [1] "code.4" "code.1" "code.3" "code.2"
# [[6]]
# [1] "code.4" "code.1" "code.2" "code.3"

从这里，您可以执行以下操作之一：

### 1, rename-in-copy
for (ord in manyorders) {
  tmpdat <- `colnames<-`(dat2, ord) # copies and renames in one line ... code-golf
  # ... your code
}

### 2, rename in place
for (ord in manyorders) {
  colnames(dat2) <- ord
  # ... your code
}

### 3, lapply, effectively rename-in-copy
all_results <- lapply(manyorders, function(ord) {
  tmpdat <- `colnames<-`(dat2, ord) # copies and renames in one line ... code-golf
  # ... your code, ending in ...
  fit <- eBayes(fit)
  y1 <- topTable(fit, coef="TreatT", n=nrow(genes), adjust.method="BH", genelist=genes)
  list(fit = fit, y1 = y1)
})

最后一个允许您查看任何运行的fit 和y1 组件，以高效的方式生成。