r - 如何使用测试统计分布对列名进行 1000 次排列?
问题描述
假设我有一个这样的矩阵
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 个列名排列的分布。
请指教
解决方案
列名的随机排序很简单:
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
组件,以有效的方式生成它。
推荐阅读
- python - 如何迭代地将编号对象添加到 python 列表中?
- java - Intellij Idea社区不识别JUnit Jupiter的注解@DisplayName
- eclipse - Eclipse 2020-03 (4.15.0) 中的 JBoss Tools 编辑器损坏
- r - R 包:如何使用 Source() 链接函数
- reporting-services - 不同高度时的并排表格布局问题
- sql - 如何在另一个加入新字段的结果上使用加入
- c++ - 为什么我的函数在重复调用时会崩溃?
- java - 如何使用 Java 函数编写程序来确定用户输入的 4 个整数是形成正方形还是矩形?
- windows - 如何从 docker 中提取 selenium/hub?
- jenkins - groovy 模块可以访问基于 jenkinsfile 的部署管道中的全局变量,这样做是个好主意吗?