【问题标题】:How can you loop this higher-order function in R?你如何在 R 中循环这个高阶函数?
【发布时间】:2015-03-20 04:00:56
【问题描述】:

这个问题与我收到的here 的回复有关,它带有来自thelatemail 的一个不错的小功能。 我使用的数据框不是最佳的,但这是我所拥有的,我只是试图在所有行中循环这个函数。

这是我的 df

dput(SO_Example_v1)
structure(list(Type = structure(c(3L, 1L, 2L), .Label = c("Community", 
"Contaminant", "Healthcare"), class = "factor"), hosp1_WoundAssocType = c(464L, 
285L, 24L), hosp1_BloodAssocType = c(73L, 40L, 26L), hosp1_UrineAssocType = c(75L, 
37L, 18L), hosp1_RespAssocType = c(137L, 77L, 2L), hosp1_CathAssocType = c(80L, 
34L, 24L), hosp2_WoundAssocType = c(171L, 115L, 17L), hosp2_BloodAssocType = c(127L, 
62L, 12L), hosp2_UrineAssocType = c(50L, 29L, 6L), hosp2_RespAssocType = c(135L, 
142L, 6L), hosp2_CathAssocType = c(95L, 24L, 12L)), .Names = c("Type", 
"hosp1_WoundAssocType", "hosp1_BloodAssocType", "hosp1_UrineAssocType", 
"hosp1_RespAssocType", "hosp1_CathAssocType", "hosp2_WoundAssocType", 
"hosp2_BloodAssocType", "hosp2_UrineAssocType", "hosp2_RespAssocType", 
"hosp2_CathAssocType"), class = "data.frame", row.names = c(NA, 
-3L))
####################
#what it looks like#
####################
require(dplyr)
df <- tbl_df(SO_Example_v1)
head(df)
         Type hosp1_WoundAssocType hosp1_BloodAssocType hosp1_UrineAssocType
1  Healthcare                  464                   73                   75
2   Community                  285                   40                   37
3 Contaminant                   24                   26                   18
Variables not shown: hosp1_RespAssocType (int), hosp1_CathAssocType (int), hosp2_WoundAssocType
  (int), hosp2_BloodAssocType (int), hosp2_UrineAssocType (int), hosp2_RespAssocType (int),
  hosp2_CathAssocType (int)

我的功能是对df$Type 中的所有类别执行chisq.test。理想情况下,如果单元格数小于 5,该函数应该切换到 fisher.test(),但这是一个单独的问题(不过,想出如何做到这一点的人会获得额外的布朗尼积分)。

这是我用来逐行执行的函数

func <- Map(
  function(x,y) {
    out <- cbind(x,y)
    final <- rbind(out[1,],colSums(out[2:3,]))
    chisq <- chisq.test(final,correct=FALSE)
    chisq$p.value
  },
  SO_Example_v1[grepl("^hosp1",names(SO_Example_v1))],
  SO_Example_v1[grepl("^hosp2",names(SO_Example_v1))] 
)
func

但理想情况下,我希望它是这样的

for(i in 1:nrow(df)){func}

但这不起作用。另一个钩子是,例如,当第二行被占用时,final 调用看起来像这样

func <- Map(
  function(x,y) {
    out <- cbind(x,y)
    final <- rbind(out[2,],colSums(out[c(1,3),]))
    chisq <- chisq.test(final,correct=FALSE)
    chisq$p.value
  },
  SO_Example_v1[grepl("^hosp1",names(SO_Example_v1))],
  SO_Example_v1[grepl("^hosp2",names(SO_Example_v1))] 
)
func

所以函数应该理解它为out[x,] 占用的单元格计数必须从colSums() 中排除。这个data.frame 只有 3 行,所以很简单,但我尝试将此函数应用于我拥有的单独的 data.frame,它包含 >200 行,所以能够以某种方式循环它会很好。

任何帮助表示赞赏。

干杯

【问题讨论】:

  • 为什么不对整个数据集执行卡方列联表?看起来您有一个 3 个类别 x 11 个类别的列联表。
  • 谢谢。原因是因为我想要用于子组比较的信息并在 3x11 列联表上运行它只会让我知道组之间存在差异,而不是哪些不同。当然,我可以查看百分比以了解正在发生的事情,但我对子组比较特别感兴趣。

标签: r function loops chi-squared


【解决方案1】:

你错过了两件事:

  1. 选择行 i 并选择除此行之外的所有要使用的行 u[i]u[-i]
  2. 如果一个项目的长度与给 Map 的其他项目的长度不同,它将被回收,这是该语言的一个非常普遍的属性。然后,您只需向函数添加一个参数,该参数对应于您要反对其他行的行,它将为传递的向量的所有项目回收。

以下内容满足您的要求

    # the function doing the stats
    FisherOrChisq <- function(x,y,lineComp) {
        out <- cbind(x,y)
        final <- rbind(out[lineComp,],colSums(out[-lineComp,]))
        test <- chisq.test(final,correct=FALSE)

        return(test$p.value)
    }

    # test of the stat function
    FisherOrChisq(SO_Example_v1[grep("^hosp1",names(SO_Example_v1))[1]],
    SO_Example_v1[grep("^hosp2",names(SO_Example_v1))[1]],2)

    # making the loop
    result <- c()
    for(type in SO_Example_v1$Type){
        line <- which(SO_Example_v1$Type==type)
        res <- Map(FisherOrChisq,
                    SO_Example_v1[grepl("^hosp1",names(SO_Example_v1))],
                    SO_Example_v1[grepl("^hosp2",names(SO_Example_v1))],
                    line
                )
        result <- rbind(result,res)
    }
    colnames(result) <- gsub("^hosp[0-9]+","",colnames(result))
    rownames(result) <- SO_Example_v1$Type

也就是说,您所做的是非常繁重的多重测试。我会非常谨慎地使用相应的 p 值,您至少需要使用多重测试校正,例如建议的 here

【讨论】:

  • 那里的编程很不错,结果也很漂亮。谢谢你。但是,如果您查看我提供的第一个函数的结果,计算出的 p 值与您的结果不同。例如 WoundAssocType 的正确列联表如下:[1,] 464 171 [2,] 309 132。并且运行chisq.test(cbind(c(464, 309), c(171, 132)), correct=F),也会为您的计算返回不同的 p 值。知道为什么吗?回复:多次测试,我知道这个问题,但不完全同意,但没有空间讨论,我们可以在以后的帖子中讨论。
  • 对于医疗保健与其他人的总和,(您提供的列联表),我的代码使用相同的列联表并返回 p 值 0.3014186,因为它当然使用了 Fisher 检验,小于5分是你的分界线吗?这里我们有 2...
  • chisq.test(cbind(c(464, 309), c(171, 132)), correct=F) 返回 pvalue 0.2815。截止值是列联表中每个单元格的计数。如果一个单元格 then fisher.test()。当所有单元格都大于 100 时,为什么您的代码会分配 fisher.test()?谢谢你的一切。
  • 已更正。如果列联表中的任何单元格小于 5,则使用 fisher.test
  • 我会为你提供赏金,因为这几乎是我所要求的(顺便谢谢你!),但我注意到,该函数将运行 fisher.test( ) 在一行的所有列联表上,如果任何一个列联表的单元格计数
猜你喜欢
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 2015-09-19
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
相关资源
最近更新 更多