4

在这里，我有下面的数据集。在 SAS 系统中，要在我所做的任何事情中获得分组结果，例如均值、GLM 或 REG，我可以将其设为：

proc sort; by B;
proc glm;
class A;
model C=A; by B;
run;

然后，我可以在 B 级以内或按 B 级获得 GLM 结果。但是，我不知道如何在 R 系统中使用 like'by' 进行分组。您可能想建议我使用>subset()，但是，例如，当您有 10 级 B 时，这条线将非常困难。您可能想学习我，但只需要 anova 分析，还需要回归和平均。这里有人可以帮助我吗？

raw data

A   B   C
a   a   0.47
a   b   0.88
a   c   2.32
a   d   3.26
a   a   0.93
a   b   1.86
a   c   3.22
a   d   0.92
a   a   0.45
a   b   0.92
a   c   2.31
a   d   3.24
b   a   0.91
b   b   1.84
b   c   3.27
b   d   0.86
b   a   0.47
b   b   0.90
b   c   2.33
b   d   3.19
b   a   0.92
b   b   1.84
b   c   3.25
b   d   0.93
c   a   0.45
c   b   0.92
c   c   2.33
c   d   3.08
c   a   0.93
c   b   1.86
c   c   3.25
c   d   0.93
c   a   0.47
c   b   0.90
c   c   2.26
c   d   3.09
4

3 回答 3

2

基于@Brian Diggs 的回答，但使用 R 基函数。我还使用了 Brian 的数据集

models <- lapply(split(dat, dat$B), function(x) glm(C~A, data=x)) 
do.call(rbind, lapply(models, function(y) y$coef))
  (Intercept)         Ab            Ac
a   0.6166667  0.1500000  1.802585e-17
b   1.2200000  0.3066667  6.666667e-03
c   2.6166667  0.3333333 -3.333333e-03
d   2.4733333 -0.8133333 -1.066667e-01

编辑 1：P.values 的结果（备选方案 1）

models <- lapply(split(dat, dat$B), function(x) glm(C~A, data=x)) #the same as above 
Coef <- lapply(models, function(y) y$coef) # the same as above
Pval <-  lapply(models, function(z) summary(z)$coefficients[, 'Pr(>|t|)'])
Result <- cbind(do.call(rbind, Coef), do.call(rbind, Pval))
colnames(Result)[4:6] <- paste('P.val', colnames(Result)[4:6])
Result

 (Intercept)         Ab            Ac P.val (Intercept)  P.val Ab  P.val Ac
a   0.6166667  0.1500000  1.802585e-17       0.007088207 0.5167724 1.0000000
b   1.2200000  0.3066667  6.666667e-03       0.008446002 0.5191737 0.9886090
c   2.6166667  0.3333333 -3.333333e-03       0.000151772 0.4762936 0.9941859
d   2.4733333 -0.8133333 -1.066667e-01       0.016803317 0.4744404 0.9235623

编辑 2：具有 P.values 的结果（备选方案 2）

models <- lapply(split(dat, dat$B), function(x) glm(C~A, data=x)) #the same as above
do.call(rbind, lapply(models, function(z) summary(z)$coefficients))
                 Estimate Std. Error       t value    Pr(>|t|)
(Intercept)  6.166667e-01  0.1540202  4.003804e+00 0.007088207
Ab           1.500000e-01  0.2178175  6.886500e-01 0.516772358
Ac           1.802585e-17  0.2178175  8.275668e-17 1.000000000
(Intercept)  1.220000e+00  0.3167661  3.851423e+00 0.008446002
Ab           3.066667e-01  0.4479749  6.845622e-01 0.519173660
Ac           6.666667e-03  0.4479749  1.488179e-02 0.988608995
(Intercept)  2.616667e+00  0.3103164  8.432253e+00 0.000151772
Ab           3.333333e-01  0.4388537  7.595545e-01 0.476293582
Ac          -3.333333e-03  0.4388537 -7.595545e-03 0.994185937
(Intercept)  2.473333e+00  0.7538543  3.280917e+00 0.016803317
Ab          -8.133333e-01  1.0661110 -7.628974e-01 0.474440418
Ac          -1.066667e-01  1.0661110 -1.000521e-01 0.923562285
于 2012-10-12T18:04:22.463 回答
2

使您的数据集更容易导入到 R：

dat <-
structure(list(A = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L), .Label = c("a", 
"b", "c"), class = "factor"), B = structure(c(1L, 2L, 3L, 4L, 
1L, 2L, 3L, 4L, 1L, 2L, 3L, 4L, 1L, 2L, 3L, 4L, 1L, 2L, 3L, 4L, 
1L, 2L, 3L, 4L, 1L, 2L, 3L, 4L, 1L, 2L, 3L, 4L, 1L, 2L, 3L, 4L
), .Label = c("a", "b", "c", "d"), class = "factor"), C = c(0.47, 
0.88, 2.32, 3.26, 0.93, 1.86, 3.22, 0.92, 0.45, 0.92, 2.31, 3.24, 
0.91, 1.84, 3.27, 0.86, 0.47, 0.9, 2.33, 3.19, 0.92, 1.84, 3.25, 
0.93, 0.45, 0.92, 2.33, 3.08, 0.93, 1.86, 3.25, 0.93, 0.47, 0.9, 
2.26, 3.09)), .Names = c("A", "B", "C"), class = "data.frame", row.names = c(NA, 
-36L))

大部分 SAS 通过组处理映射到拆分-应用-组合方法（将数据分成几部分，对每个部分做一些事情，以某种方式将这些部分重新组合在一起）。在这种情况下，模型的结果是对象（列表），“组合”多个模型的自然方法是将它们放在列表中。

library("plyr")
models <- dlply(dat, .(B), function(DF) glm(C~A, data=DF))

models现在是一个列表，其中的每个元素都是 的glm子集上 a的结果dim。

> models
$a

Call:  glm(formula = C ~ A, data = DF)

Coefficients:
(Intercept)           Ab           Ac  
  6.167e-01    1.500e-01    6.799e-17  

Degrees of Freedom: 8 Total (i.e. Null);  6 Residual
Null Deviance:      0.472 
Residual Deviance: 0.427        AIC: 6.107 

$b

Call:  glm(formula = C ~ A, data = DF)

Coefficients:
(Intercept)           Ab           Ac  
   1.220000     0.306667     0.006667  

Degrees of Freedom: 8 Total (i.e. Null);  6 Residual
Null Deviance:      1.99 
Residual Deviance: 1.806        AIC: 19.09 

$c

Call:  glm(formula = C ~ A, data = DF)

Coefficients:
(Intercept)           Ab           Ac  
   2.616667     0.333333    -0.003333  

Degrees of Freedom: 8 Total (i.e. Null);  6 Residual
Null Deviance:      1.958 
Residual Deviance: 1.733        AIC: 18.72 

$d

Call:  glm(formula = C ~ A, data = DF)

Coefficients:
(Intercept)           Ab           Ac  
     2.4733      -0.8133      -0.1067  

Degrees of Freedom: 8 Total (i.e. Null);  6 Residual
Null Deviance:      11.4 
Residual Deviance: 10.23        AIC: 34.69 

attr(,"split_type")
[1] "data.frame"
attr(,"split_labels")
  B
1 a
2 b
3 c
4 d

一次从所有模型中提取信息遵循相同的范例：

> ldply(models, coefficients)
  B (Intercept)         Ab            Ac
1 a   0.6166667  0.1500000  6.798700e-17
2 b   1.2200000  0.3066667  6.666667e-03
3 c   2.6166667  0.3333333 -3.333333e-03
4 d   2.4733333 -0.8133333 -1.066667e-01
于 2012-10-12T17:37:47.613 回答
0

也许聚合函数就是你要找的

aggregate(data$C, by=list(data$A), FUN=sum)

这将按第一列对您的数据进行分组，并将 C 列折叠为第一列中每个组的总和

于 2012-10-12T16:50:29.357 回答