r-使用方差分析结果列表中的p值作为条件

Question

我想运行一系列kruskal.test，然后运行一个dunn_test，其中kruskal是重要的。然后打印重要dunn_tests的结果。

library(rstatix)
library(purrr)
iris<- iris
for(i in 1:4){
a<- colnames(iris)
anova<- map(names(iris)[1:4], ~ kruskal.test(reformulate('Species', response=.x), data=iris))
post<- dunn_test(iris, as.formula(paste(a[i], a[5], sep="~")), p.adjust.method= "bonferroni")

if(any(post$p.adj < 0.05 & anova[[i]]$p.value < 0.05)){
print(a)
print(setNames(post, a))
}
}

我想知道，当同一变量的kruskal测试是显著的时，我怎么能只打印显著的Dunnstest。

我意识到这不会给我想要的结果，因为帖子和anova项目是两个独立的对象，这只是我最初的想法。使用anova[[i]]$p.value还会导致下标越界错误。

Answer 1

这里有一种使用函数来获得所需的临时结果的方法。它的工作原理类似于循环。我知道你对as.formula很熟悉，我认为这是最好的方法。我在Kruskal-Wallis测试和Dunn测试中使用了相同的公式。该函数只使用要测试的列名，并根据该列名创建公式。如果您想返回所有的事后测试(如果任何一个比较是显著的(，any只需进行一个小的更改。输出want捕获具有至少一个显著成对比较的所有事后测试。

除了代码之外，创建与R中的函数同名的对象名称并不是一个好的做法。anova是R中的基函数，因此您不希望将结果命名为anova。您还可以在循环中一次又一次地创建对象a。这是同一项任务，所以它应该在您的循环之外。这不是什么大不了的事，但如果你开始使用更大的数据集，它可能会累积起来。

library(rstatix)
#> 
#> Attaching package: 'rstatix'
#> The following object is masked from 'package:stats':
#> 
#>     filter
iris <- iris
a <- colnames(iris)
a <- a[1:4]
myfunc <- function(y, x){
form <- as.formula(paste0(y, " ~ ", x))
res <- rstatix::kruskal_test(form, data = iris)

if (res$p < 0.05){
post <- rstatix::dunn_test(data = iris, form, p.adjust.method= "bonferroni")

if(any(post$p.adj < 0.05)) return(post)
}
}
want <- lapply(a, FUN = myfunc, x = "Species")
want <- do.call(rbind, want)

如果你只想要有意义的特定成对比较，因为dunn_test的输出是data.frame，你可以用小于0.05的调整后的p值和小于0.05的return来子集。

myfunc <- function(y, x){
form <- as.formula(paste0(y, " ~ ", x))
res <- rstatix::kruskal_test(form, data = iris)

if (res$p < 0.05){
post <- rstatix::dunn_test(data = iris, form, p.adjust.method= "bonferroni")
post <- post[post$p.adj < 0.05 ,]
return(post)
}
}

作为一个循环，这将是：

res <- list()
post <- list()
post_sig <- list()
for(i in 1:4){
form <- as.formula(paste0(a[i], " ~ Species"))
res[[i]] <- rstatix::kruskal_test(form, data = iris)
post[[i]] <- rstatix::dunn_test(data = iris, form, p.adjust.method= "bonferroni")

if (res[[i]]$p < 0.05 & any(post[[i]]$p.adj < 0.05)){
print(post[[i]])
## could be blanks in the list here as the number i is used
## e.g. if var 2 was not significant, post_sig[[2]] is blank
post_sig[[i]] <- post[[i]]
}
}
## if you want results as a dataframe
res <- do.call(rbind, res)
post <- do.call(rbind, post)
post_sig <- do.call(rbind, post_sig)

Answer 2

据我所知，您的代码大部分都能工作。唯一的主要问题是您滥用了any()函数。你把第二个括号放错地方了。修复并重新排列代码，我运行了这个：

library(rstatix)
library(purrr)
iris<- iris
a<- colnames(iris)
anova<- map(names(iris)[1:4], ~ kruskal.test(reformulate('Species', response=.x), data=iris))
for(i in 1:4){
post<- dunn_test(iris, as.formula(paste(a[i], a[5], sep="~")), p.adjust.method= "bonferroni")

if(any(post$p.adj < 0.05) & anova[[i]]$p.value < 0.05){
print(a[i])
print(post)
}
}

只要kruskal测试的相关p值小于0.05，并且dunn_test中至少有一个比较有效，这段代码就会一个接一个地打印出被测试变量的名称和dunn_test的输出。