通过 anova 与 purrr 或 dplyr 比较模型答案

【问题标题】：Compare models via anova with purrr or dplyr通过 anova 与 purrr 或 dplyr 比较模型
【发布时间】：2018-11-15 01:29:43
【问题描述】：

我想拟合几个模型并通过anova 函数进行比较。有一个 link 描述了如何使用 do 完成此操作，但 anova 调用对我不起作用。请注意，我现在使用 purrr::map 但它不会改变我得到的错误（见下文）

这是我的代码：

out = df %>% 
  group_by(condition,File) %>%
  nest() %>% 
  mutate(fit_nls = purrr::map(data, ~ minpack.lm::nlsLM(rating ~ a*cd^b + d,data=.)),
         fit_lm = purrr::map(data, ~ lm(rating ~ cd,data=.)))

out %>% do(aov = stats::anova(.$fit_nls, .$fit_lm))

Error in UseMethod("anova") : 
  no applicable method for 'anova' applied to an object of class "list"*

会话信息：

R version 3.4.2 (2017-09-28)
Platform: x86_64-apple-darwin15.6.0 (64-bit)
Running under: macOS High Sierra 10.13.4

Matrix products: default
BLAS: /System/Library/Frameworks/Accelerate.framework/Versions/A/Frameworks/vecLib.framework/Versions/A/libBLAS.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/3.4/Resources/lib/libRlapack.dylib

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] readbulk_1.1.0   assertr_2.5      afex_0.20-2      emmeans_1.2      lme4_1.1-17      Matrix_1.2-11    naniar_0.2.0     car_3.0-0       
 [9] knitr_1.20       kableExtra_0.8.0 broom_0.4.4      bindrcpp_0.2.2   effects_4.0-1    carData_3.0-1    QuantPsyc_1.5    MASS_7.3-47     
[17] boot_1.3-20      jmv_0.8.6.2      forcats_0.3.0    stringr_1.3.0    dplyr_0.7.4      purrr_0.2.4      readr_1.1.1      tidyr_0.8.0     
[25] tibble_1.4.2     ggplot2_2.2.1    tidyverse_1.2.1  pacman_0.4.6    

loaded via a namespace (and not attached):
  [1] readxl_1.1.0         backports_1.1.2      plyr_1.8.4           igraph_1.2.1         lazyeval_0.2.1       splines_3.4.2       
  [7] crosstalk_1.0.0      TH.data_1.0-8        rstantools_1.5.0     inline_0.3.14        digest_0.6.15        htmltools_0.3.6     
 [13] rsconnect_0.8.8      lmerTest_3.0-1       magrittr_1.5         openxlsx_4.0.17      brms_2.2.0           modelr_0.1.2        
 [19] matrixStats_0.53.1   xts_0.10-2           sandwich_2.4-0       colorspace_1.3-2     rvest_0.3.2          haven_1.1.1         
 [25] crayon_1.3.4         jsonlite_1.5         bindr_0.1.1          survival_2.41-3      zoo_1.8-1            glue_1.2.0          
 [31] gtable_0.2.0         rstan_2.17.3         abind_1.4-5          scales_0.5.0         mvtnorm_1.0-7        miniUI_0.1.1        
 [37] Rcpp_0.12.16         viridisLite_0.3.0    xtable_1.8-2         foreign_0.8-69       stats4_3.4.2         StanHeaders_2.17.2  
 [43] survey_3.33-2        DT_0.4               htmlwidgets_1.2      httr_1.3.1           threejs_0.3.1        modeltools_0.2-21   
 [49] pkgconfig_2.0.1      loo_2.0.0            nnet_7.3-12          utf8_1.1.3           tidyselect_0.2.4     labeling_0.3        
 [55] rlang_0.2.0          reshape2_1.4.3       later_0.7.1          munsell_0.4.3        cellranger_1.1.0     tools_3.4.2         
 [61] cli_1.0.0            jmvcore_0.8.5        ggridges_0.5.0       evaluate_0.10.1      yaml_2.1.18          coin_1.2-2          
 [67] visdat_0.1.0         nlme_3.1-131         mime_0.5             xml2_1.2.0           compiler_3.4.2       bayesplot_1.5.0     
 [73] shinythemes_1.1.1    rstudioapi_0.7       curl_3.2             stringi_1.1.6        highr_0.6            Brobdingnag_1.2-5   
 [79] lattice_0.20-35      psych_1.8.4          nloptr_1.0.4         markdown_0.8         shinyjs_1.0          pillar_1.2.2        
 [85] bridgesampling_0.4-0 estimability_1.3     data.table_1.10.4-3  httpuv_1.4.0         R6_2.2.2             promises_1.0.1      
 [91] gridExtra_2.3        rio_0.5.10           codetools_0.2-15     colourpicker_1.0     gtools_3.5.0         assertthat_0.2.0    
 [97] rprojroot_1.3-2      rjson_0.2.15         minpack.lm_1.2-1     shinystan_2.5.0      mnormt_1.5-5         multcomp_1.4-8      
[103] parallel_3.4.2       hms_0.4.2            grid_3.4.2           coda_0.19-1          minqa_1.2.4          rmarkdown_1.9       
[109] numDeriv_2016.8-1    shiny_1.0.5          lubridate_1.7.4      base64enc_0.1-3      dygraphs_1.1.1.4

这是数据：

dput(df)
structure(list(condition = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L), .Label = c("new", "standard"
), class = "factor"), File = c("001", "001", "001", "001", "001", 
"001", "001", "001", "002", "002", "002", "002", "002", "002", 
"002", "002", "003", "003", "003", "003", "003", "003", "003", 
"003", "004", "004", "004", "004", "004", "004", "004", "004", 
"005", "005", "005", "005", "005", "005", "005", "005", "006", 
"006", "006", "006", "006", "006", "006", "006", "007", "007", 
"007", "007", "007", "007", "007", "007", "008", "008", "008", 
"008", "008", "008", "008", "008", "009", "009", "009", "009", 
"009", "009", "009", "009", "010", "010", "010", "010", "010", 
"010", "010", "010", "011", "011", "011", "011", "011", "011", 
"011", "011", "012", "012", "012", "012", "012", "012", "012", 
"012", "013", "013", "013", "013", "013", "013", "013", "013", 
"014", "014", "014", "014", "014", "014", "014", "014", "015", 
"015", "015", "015", "015", "015", "015", "015", "016", "016", 
"016", "016", "016", "016", "016", "016", "017", "017", "017", 
"017", "017", "017", "017", "017", "019", "019", "019", "019", 
"019", "019", "019", "019", "020", "020", "020", "020", "020", 
"020", "020", "020", "021", "021", "021", "021", "021", "021", 
"021", "021", "022", "022", "022", "022", "022", "022", "022", 
"022", "023", "023", "023", "023", "023", "023", "023", "023", 
"024", "024", "024", "024", "024", "024", "024", "024", "025", 
"025", "025", "025", "025", "025", "025", "025", "026", "026", 
"026", "026", "026", "026", "026", "026", "027", "027", "027", 
"027", "027", "027", "027", "027", "028", "028", "028", "028", 
"028", "028", "028", "028", "029", "029", "029", "029", "029", 
"029", "029", "029", "031", "031", "031", "031", "031", "031", 
"031", "031", "032", "032", "032", "032", "032", "032", "032", 
"032", "033", "033", "033", "033", "033", "033", "033", "033", 
"034", "034", "034", "034", "034", "034", "034", "034", "035", 
"035", "035", "035", "035", "035", "035", "035", "036", "036", 
"036", "036", "036", "036", "036", "036"), rating = c(15, 12, 
7.69230769230769, 11, 6.66666666666667, 5, 8.33333333333333, 
16, 5, 8, 11, 20, 6, 13, 18, 5, 14, 7.69230769230769, 11, 6.25, 
14, 7.14285714285714, 5.55555555555556, 18, 7, 7, 12, 19, 2, 
3, 10, 20, 20, 11, 6.66666666666667, 5, 25, 12, 6.66666666666667, 
4, 50, 30, 5, 11, 7.5, 2, 1, 15, 6.66666666666667, 8.33333333333333, 
20, 3.33333333333333, 5, 50, 15, 11, 9.25, 7, 10.5, 11, 18, 5.5, 
13, 2, 4, 12, 16, 6.66666666666667, 25, 5, 8.33333333333333, 
17, 6, 15, 8, 4, 13, 9, 30, 11, 12, 90, 2, 3.33333333333333, 
5, 20, 6.25, 30, 12, 30, 7, 4, 4, 6, 20, 15, 13, 12, 6.66666666666667, 
25, 6.66666666666667, 18, 5, 9.09090909090909, 8, 20, 14, 2, 
11, 6, 9, 12, 6.66666666666667, 45, 18, 2, 5, 12, 25, 4, 20, 
1.5, 13, 18, 9, 9, 7, 2, 30, 15, 5.55555555555556, 5, 11, 7.69230769230769, 
20, 3.33333333333333, 14, 7.14285714285714, 7.69230769230769, 
5.88235294117647, 11, 16, 18, 8.33333333333333, 7, 8, 3, 11, 
4, 14, 20, 21, 7.14285714285714, 14, 20, 12, 6.66666666666667, 
8.33333333333333, 7.69230769230769, 11, 6, 12, 18, 4, 7, 13, 
2, 11, 4, 13, 6.66666666666667, 30, 5, 30, 20, 2.85714285714286, 
20, 11, 17, 2, 6, 13, 9, 4, 2.5, 7.14285714285714, 35, 5.26315789473684, 
4, 11, 16, 25, 5, 14, 11, 2, 8, 7, 18, 20, 16, 20, 8.33333333333333, 
4.76190476190476, 11, 7.14285714285714, 5.88235294117647, 12, 
17, 8, 23, 5, 10.5, 6, 11, 2.6, 6.66666666666667, 12, 6.25, 30, 
15, 3.33333333333333, 6.66666666666667, 25, 5.55555555555556, 
7.69230769230769, 4.34782608695652, 15, 20, 11, 8.33333333333333, 
12, 2, 17, 7, 15, 5, 20, 11, 6, 2.5, 3.33333333333333, 30, 15, 
11, 6.66666666666667, 50, 4, 7, 26, 4, 10, 3, 13, 20, 9, 11, 
6.66666666666667, 2.5, 40, 13, 8.33333333333333, 5, 18, 13, 7, 
1, 20, 11, 1, 7, 17), cd = c(57.2, 32, 3.2, 17.9, 1.8, 1, 5.7, 
100, 1, 5.7, 17.9, 100, 3.2, 32, 57.2, 1.8, 32, 5.7, 17.9, 3.2, 
57.2, 1.8, 1, 100, 3.2, 5.7, 32, 57.2, 1, 1.8, 17.9, 100, 57.2, 
17.9, 1.8, 3.2, 100, 32, 5.7, 1, 100, 57.2, 3.2, 17.9, 5.7, 1.8, 
1, 32, 3.2, 5.7, 57.2, 1, 1.8, 100, 32, 17.9, 5.7, 3.2, 17.9, 
32, 100, 1.8, 57.2, 1, 1, 17.9, 32, 3.2, 100, 1.8, 5.7, 57.2, 
1.8, 57.2, 3.2, 1, 32, 5.7, 100, 17.9, 17.9, 100, 1, 1.8, 3.2, 
32, 5.7, 57.2, 17.9, 100, 5.7, 1.8, 1, 3.2, 57.2, 32, 32, 17.9, 
1.8, 100, 3.2, 57.2, 1, 5.7, 3.2, 100, 57.2, 1, 17.9, 1.8, 5.7, 
32, 5.7, 100, 32, 1, 1.8, 17.9, 57.2, 3.2, 100, 1, 32, 57.2, 
5.7, 17.9, 3.2, 1.8, 100, 32, 3.2, 1.8, 17.9, 5.7, 57.2, 1, 32, 
1.8, 3.2, 1, 17.9, 57.2, 100, 5.7, 3.2, 5.7, 1.8, 17.9, 1, 32, 
57.2, 100, 1.8, 57.2, 100, 32, 1, 5.7, 3.2, 17.9, 3.2, 32, 100, 
1.8, 5.7, 57.2, 1, 17.9, 1.8, 17.9, 5.7, 100, 3.2, 57.2, 32, 
1, 100, 17.9, 57.2, 1, 3.2, 32, 5.7, 1.8, 1, 5.7, 100, 3.2, 1.8, 
17.9, 32, 57.2, 1.8, 32, 17.9, 1, 5.7, 3.2, 57.2, 100, 57.2, 
100, 5.7, 1, 17.9, 3.2, 1.8, 32, 57.2, 5.7, 100, 1.8, 17.9, 3.2, 
32, 1, 5.7, 17.9, 3.2, 100, 32, 1, 1.8, 57.2, 1.8, 3.2, 1, 57.2, 
100, 17.9, 5.7, 32, 1, 57.2, 5.7, 32, 1.8, 100, 17.9, 3.2, 1, 
1.8, 57.2, 32, 17.9, 5.7, 100, 3.2, 3.2, 100, 1.8, 17.9, 1, 32, 
57.2, 5.7, 17.9, 3.2, 1, 100, 32, 5.7, 1.8, 57.2, 32, 5.7, 1.8, 
100, 17.9, 1, 3.2, 57.2)), .Names = c("condition", "File", "rating", 
"cd"), row.names = c(NA, -272L), class = "data.frame")

【问题讨论】：

标签： r dplyr anova purrr

【解决方案1】：

您有两个模型列表；如果要比较每组中的两个模型，则需要循环遍历两个列表并调用anova；一种方法是使用来自purrr 包的map2：

out %>% 
    mutate(aov = map2(fit_nls, fit_lm, stats::anova)) %>% 
    pull(aov) %>% head(2)

#[[1]]
#Analysis of Variance Table

#Model 1: rating ~ a * cd^b + d
#Model 2: rating ~ cd
#  Res.Df Res.Sum Sq Df  Sum Sq F value    Pr(>F)    
#1      5     1.2729                                 
#2      6    17.8237 -1 -16.551  65.014 0.0004751 ***
#---
#Signif. codes:  
#0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

#[[2]]
#Analysis of Variance Table

#Model 1: rating ~ a * cd^b + d
#Model 2: rating ~ cd
#  Res.Df Res.Sum Sq Df  Sum Sq F value   Pr(>F)   
#1      5      2.838                               
#2      6     23.496 -1 -20.658  36.395 0.001802 **
#---
#Signif. codes:  
#0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

【讨论】：

谢谢！我仍然想知道为什么链接中的代码似乎可以正常工作，尽管我做了同样的事情（除了 purrr:map 而不是 do）。
链接起作用的原因是models 被分组，所以链接do 之后将anova 应用于每个组。在您的情况下，nest ungroup 数据框，do 然后将 anova 应用于整个数据框而不是按组。

【解决方案2】：

您正在混合两种相互竞争的方法：（即将）弃用的group_by + do 和（现在）推荐的nest + mutate - map。您应该使用@Psidom 的解决方案，但是如果您想了解为什么您的代码不起作用，那是因为out 不再“分组”，而且.$fit_nls 是一个列表，您需要采取每个组的第一个（也是唯一一个）元素：

out %>%
  group_by(condition, File) %>%
  do(aov = stats::anova(.$fit_nls[[1]], .$fit_lm[[1]]))

（编辑：在他们自己的答案下没有看到 Psidom 的评论，这几乎可以概括）。

【讨论】：