0

我有一个glm模型和一个汇总数据集,需要我将coefficients,standard errorp.value从模型的摘要绑定到汇总数据集。例如,我使用了mtcars数据集。我添加columns到最终的联合数据集中以模拟我想要coefficients, standard errors, and p-values放置的位置。对于模型中未显示的基本值,我想1在. 我怎么能做到这一切?coefficientsintercept, standard errors and p-value

library(tidyverse)

mtcars <- as_tibble(mtcars)

mtcars$cyl <- as.factor(mtcars$cyl)
mtcars$gear <- as.factor(mtcars$gear)

#run model
model1 <- glm(mpg ~ cyl + gear, data = mtcars)
summary(model1)

#start developing summarized data set
mtcars_wght <- mtcars %>%
  group_by(cyl) %>%
  rename(level = cyl) %>%
  summarise("sum_weight" = sum(wt)) %>%
  mutate("variable" = "cyl") 

mtcars_gear <- mtcars %>%
  group_by(gear) %>%
  summarise("sum_weight" = sum(wt)) %>%
  mutate("variable" = "gear") %>%
  rename(level = gear)

#make summarized data set example
mtcars_sum <- mtcars_wght %>%
  bind_rows(mtcars_gear) %>%
  mutate("coefficient" = "x", "std.error" = "y", "p_value" = "z")
4

1 回答 1

1

你在这个输出之后吗?

# A tibble: 6 x 6
  level sum_weight variable coefficient std.error    p_value
  <chr>      <dbl> <chr>          <dbl>     <dbl>      <dbl>
1 4           25.1 cyl            NA        NA    NA        
2 6           21.8 cyl            -6.66      1.63  0.000353 
3 8           56.0 cyl           -10.5       1.96  0.0000109
4 3           58.4 gear           NA        NA    NA        
5 4           31.4 gear            1.32      1.93  0.498    
6 5           13.2 gear            1.50      1.85  0.426

以下是您如何使用dplyrand来实现这一点broom

df <- rbind(mtcars_wght, mtcars_gear)
df <- df %>% mutate(
  level = paste0(variable, level)
) %>% select(-variable)

mod_summary <- model1 %>% broom::tidy()

left_join(df, mod_summary, by = c('level' = 'term')) %>% 
  mutate(variable = str_extract(level, '[a-z]+'),
         level = str_extract(level, '[0-9]+')) %>%
  rename(coefficient = estimate, p_value = p.value) %>%
  select(level, sum_weight, variable, coefficient, std.error, p_value)

编辑

如果要包含Intercept,请使用full_join而不是left_join上面。下面,我将输出保存到thesummary.

thesummary <- full_join(df, mod_summary, by = c('level' = 'term')) %>% 
  mutate(variable = str_extract(level, '[A-Za-z]+'),
         level = str_extract(level, '[0-9]+')) %>%
  rename(coefficient = estimate, p_value = p.value) %>%
  select(level, sum_weight, variable, coefficient, std.error, p_value)

要分配1缺失值,仅针对最后 4 列,请执行以下操作:

cbind(thesummary[,1:2], apply(thesummary[,3:6], 2, function(x) ifelse(is.na(x), 1, x)))

这是输出:

  level sum_weight  variable coefficient std.error   p_value
1     4      25.14       cyl           1         1         1
2     6      21.82       cyl      -6.656     1.629 3.528e-04
3     8      55.99       cyl     -10.542     1.958 1.087e-05
4     3      58.39      gear           1         1         1
5     4      31.40      gear       1.324     1.928 4.980e-01
6     5      13.16      gear       1.500     1.855 4.257e-01
7  <NA>         NA Intercept      25.428     1.881 1.554e-13

如果你想NA用`替换每个,只需执行:

thesummary[is.na(thesummary)] <- 1

这是输出。

# A tibble: 7 x 6
  level sum_weight variable  coefficient std.error  p_value
  <chr>      <dbl> <chr>           <dbl>     <dbl>    <dbl>
1 4          25.1  cyl              1.00      1.00 1.00e+ 0
2 6          21.8  cyl             -6.66      1.63 3.53e- 4
3 8          56.0  cyl            -10.5       1.96 1.09e- 5
4 3          58.4  gear             1.00      1.00 1.00e+ 0
5 4          31.4  gear             1.32      1.93 4.98e- 1
6 5          13.2  gear             1.50      1.85 4.26e- 1
7 1           1.00 Intercept       25.4       1.88 1.55e-13
于 2018-04-05T09:57:30.623 回答