4

我想从书中重现一些计算（logit 回归）。这本书给出了一个列联表和结果。

这是表格：

在此处输入图像描述

    .


    example <- matrix(c(21,22,6,51), nrow = 2, byrow = TRUE)
    #Labels:
    rownames(example) <- c("Present","Absent")
    colnames(example) <- c(">= 55", "<55")

它给了我这个：

            >= 55 <55
    Present    21  22
    Absent      6  51

但是要使用 glm() 函数，数据必须采用以下方式：

（两列，一列“Age”，一列“Present”，用0/1填充）

    age <- c(rep(c(0),27), rep(c(1),73))
    present <- c(rep(c(0),21), rep(c(1),6), rep(c(0),22), rep(c(1),51))

    data <- data.frame(present, age)

    > data
        present age
    1         0   0
    2         0   0
    3         0   0
    .         .   .
    .         .   .
    .         .   .
    100       1   1

有没有一种简单的方法可以从表/矩阵中获取这种结构？

4

5 回答 5

2

您或许可以使用此处定义countsToCases的函数。

countsToCases(as.data.frame(as.table(example))) 
#        Var1  Var2
#1    Present >= 55
#1.1  Present >= 55
#1.2  Present >= 55
#1.3  Present >= 55
#1.4  Present >= 55
#1.5  Present >= 55
# ...

如果您愿意，您可以随时将变量重新编码为数字。

于 2016-04-08T12:34:38.090 回答
2 
reshape2::melt(example)

这会给你，

     Var1  Var2 value
1 Present >= 55    21
2  Absent >= 55     6
3 Present   <55    22
4  Absent   <55    51

您可以轻松地用于glm

于 2016-04-08T13:32:18.407 回答
1

我会去：

library(data.table)
tab <- data.table(AGED = c(1, 1, 0, 0),
                  CHD = c(1, 0, 1, 0),
                  Count = c(21, 6, 22, 51))

tabExp <- tab[rep(1:.N, Count), .(AGED, CHD)]

编辑：快速解释，因为我花了一些时间才弄清楚：

在data.tableobjects 中.N存储组的行数（如果与 分组by）或仅存储整体的行数data.table，因此在此示例中：

tab[rep(1:.N, Count)]

和

tab[rep(1:4, Count)]

最后

tab[rep(1:4, c(21, 6, 22, 51)]

是等价的。

与基础 R 相同：

tab2 <- data.frame(AGED = c(1, 1, 0, 0),
                   CHD = c(1, 0, 1, 0),
                   Count = c(21, 6, 22, 51))

tabExp2 <- tab2[rep(1:nrow(tab2), tab2$Count), c("AGED", "CHD")]
于 2016-04-08T12:30:01.300 回答
1

所以，glm并不是那么死板。部分?glm读取

 For ‘binomial’ and ‘quasibinomial’ families the response can also
 be specified as a ‘factor’ (when the first level denotes failure
 and all others success) or as a two-column matrix with the columns
 giving the numbers of successes and failures.

我假设您想测试年龄对Present/Absent. 关键是用于指定响应，如 (in psueudo-code) c(success, failure)。

所以你需要像data.frame(Age= ..., Present = ..., Absent). 从您那里执行此操作的最简单方法example是转置，然后强制转换为data.frame，并添加一列：

example_t <- as.data.frame(t(example))
example_df <- data.frame(example_t, Age=factor(row.names(example_t)))

这给了你

      Present Absent   Age
>= 55      21      6 >= 55
<55        22     51   <55

然后，您可以运行 glm：

glm(cbind(Present, Absent) ~ Age, example_df, family = 'binomial')

要得到

Call:  glm(formula = cbind(Present, Absent) ~ Age, family = "binomial",
    data = example_for_glm)

Coefficients:
(Intercept)       Age<55
      1.253       -2.094

Degrees of Freedom: 1 Total (i.e. Null);  0 Residual
Null Deviance:      18.7
Residual Deviance: -1.332e-15   AIC: 11.99

附录

您也可以通过@therimalaya 的回答到达这里。但这只是第一步

as.data.frame(as.table(example))

（只会让你分道扬镳）

    Var1  Var2 Freq
1 Present >= 55   21
2  Absent >= 55    6
3 Present   <55   22
4  Absent   <55   51

但是要真正拥有成功和失败的专栏，您需要做更多的事情。你可以tidyr用来到达那里

as.data.frame(as.table(example)) %>% tidyr::spread(Var1, Freq)

example_df和我上面的类似

  Var2 Present Absent
1 >= 55      21      6
2   <55      22     51
于 2016-04-08T20:13:13.400 回答
1

下面的代码可能看起来很长，但只有group_by()anddo()指令处理扩展数据。剩下的就是以长格式更改数据并将字符变量编码为 0 和 1。我试图从您在问题中给出的确切矩阵开始。

加载数据操作包

library(tidyr)
library(dplyr)

创建数据框

像在您的示例中一样创建一个矩阵，但避免在列名中使用“>”符号

example <- matrix(c(21,22,6,51), nrow = 2, byrow = TRUE)
rownames(example) <- c("Present","Absent")
colnames(example) <- c("above55", "below55")

将矩阵转换为数据框

example <- data.frame(example) %>%
    add_rownames("chd")

或者干脆直接创建一个数据框

data.frame(chd = c("Present", "Absent"),
           above55 = c(21,6),
           below55 = c(22,51))

重塑数据

data2 <- example %>% 
    gather(age, nrow, -chd) %>%
    # Encode chd and age as 0 or 1
    mutate(chd = ifelse(chd=="Present",1,0),
           age = ifelse(age=="above55",1,0)) %>%
    group_by(chd, age) %>%
    # Expand each variable by nrow
    do(data.frame(chd = rep(.$chd,.$nrow),
                  age = rep(.$age,.$nrow)))

head(data2)
# Source: local data frame [6 x 2]
# Groups: chd, age [1]
# 
#     chd   age
#   (dbl) (dbl)
# 1     0     0
# 2     0     0
# 3     0     0
# 4     0     0
# 5     0     0
# 6     0     0
tail(data2)
# Source: local data frame [6 x 2]
# Groups: chd, age [1]
# 
#     chd   age
#   (dbl) (dbl)
# 1     1     1
# 2     1     1
# 3     1     1
# 4     1     1
# 5     1     1
# 6     1     1

table(data2)
#        age
# chd  0  1
#   0 51  6
#   1 22 21

除了我上面评论中提到的年龄编码问题外，与您的示例相同。

于 2016-04-08T13:25:18.863 回答