代码之家  ›  专栏  ›  技术社区  ›  Tito Sanz

以整洁的方式通过TIBLE运行aov测试

  •  1
  • Tito Sanz  · 技术社区  · 7 年前

    我想使用相同的因变量在数据帧上运行线性回归。类似的问题也得到了解决 here .问题是 aov 用于实现方差分析的函数不接受 x y 作为论据(据我所知)。是否有一种方法可以整齐地执行分析?到目前为止,我已经尝试过以下方法:

    library(tidyverse)
    
    iris %>% 
      as_tibble() %>% 
      select(Sepal.Length, Species) %>% 
      mutate(foo_a = as_factor(sample(c("a", "b", "c"), nrow(.), replace = T)),
             foo_b = as_factor(sample(c("d", "e", "f"), nrow(.), replace = T))) %>% 
      map(~aov(Sepal.Length ~ .x, data = .))
    

    于2019年2月12日由 reprex package (v0.2.1)

    所需的输出是三种分析: Sepal.Length Species , 萼片。长 foo_a 最后一个呢 萼片。长 foo_b .可能吗?还是我完全错了?

    0 回复  |  直到 7 年前
        1
  •  2
  •   camille    7 年前

    一种方法是将其制成一个长形数据框架,由感兴趣的自变量分组,并使用 "many models" 方法我通常更喜欢这样的东西,而不是尝试在多个专栏中进行tidyeval。这只会让我更清楚地了解正在发生的事情。

    为了节省空间,我正在与 iris_foo ,这是通过2条变异线创建的数据。将其转换为长格式将为您提供这三个列的名称键,这三个列将用作每个列中的独立变量 aov 电话。

    library(tidyverse)
    
    iris_foo %>%
      gather(key, value, -Sepal.Length)
    
    #> # A tibble: 450 x 3
    #>    Sepal.Length key     value 
    #>           <dbl> <chr>   <chr> 
    #>  1          5.1 Species setosa
    #>  2          4.9 Species setosa
    #>  3          4.7 Species setosa
    #>  4          4.6 Species setosa
    #>  5          5   Species setosa
    #>  6          5.4 Species setosa
    #>  7          4.6 Species setosa
    #>  8          5   Species setosa
    #>  9          4.4 Species setosa
    #> 10          4.9 Species setosa
    #> # … with 440 more rows
    

    从那里,筑巢 key 并创建ANOVA模型的新列表列。这将是一个列表 aov 物体。为了简化模型的恢复,您可以删除数据列。

    aov_models <- iris_foo %>%
      gather(key, value, -Sepal.Length) %>%
      group_by(key) %>%
      nest() %>%
      mutate(model = map(data, ~aov(Sepal.Length ~ value, data = .))) %>%
      select(-data)
    
    aov_models
    #> # A tibble: 3 x 2
    #>   key     model    
    #>   <chr>   <list>   
    #> 1 Species <S3: aov>
    #> 2 foo_a   <S3: aov>
    #> 3 foo_b   <S3: aov>
    

    从那里,你可以随心所欲地使用模型。它们在列表中是可访问的 aov_models$model .打印出来的,看起来就像你期望的那样。例如,第一个模型:

    aov_models$model[[1]]
    #> Call:
    #>    aov(formula = Sepal.Length ~ value, data = .)
    #> 
    #> Terms:
    #>                    value Residuals
    #> Sum of Squares  63.21213  38.95620
    #> Deg. of Freedom        2       147
    #> 
    #> Residual standard error: 0.5147894
    #> Estimated effects may be unbalanced
    

    要查看所有模型,请致电 aov_models$model %>% map(print) 。您可能还想使用 broom 功能,例如 broom::tidy broom::glance ,这取决于您需要如何展示模型。

    aov_models$model %>%
      map(broom::tidy)
    #> [[1]]
    #> # A tibble: 2 x 6
    #>   term         df sumsq meansq statistic   p.value
    #>   <chr>     <dbl> <dbl>  <dbl>     <dbl>     <dbl>
    #> 1 value         2  63.2 31.6        119.  1.67e-31
    #> 2 Residuals   147  39.0  0.265       NA  NA       
    #> 
    #> [[2]]
    #> # A tibble: 2 x 6
    #>   term         df   sumsq meansq statistic p.value
    #>   <chr>     <dbl>   <dbl>  <dbl>     <dbl>   <dbl>
    #> 1 value         2   0.281  0.141     0.203   0.817
    #> 2 Residuals   147 102.     0.693    NA      NA    
    #> 
    #> [[3]]
    #> # A tibble: 2 x 6
    #>   term         df   sumsq meansq statistic p.value
    #>   <chr>     <dbl>   <dbl>  <dbl>     <dbl>   <dbl>
    #> 1 value         2   0.756  0.378     0.548   0.579
    #> 2 Residuals   147 101.     0.690    NA      NA
    

    或者将所有模型整理成一个数据帧,这样可以保持 钥匙 列中,您可以执行以下操作:

    aov_models %>%
      mutate(model_tidy = map(model, broom::tidy)) %>%
      unnest(model_tidy)