如何用R实现SAS t检验输出结果

1. data aaa;
   
2. 
   
3. input x g;
   
4. 
   
5. cards;
   
6. 
   
7. 12  1
   
8. 
   
9. 13  1
   
10. 
    
11. 14  1
    
12. 
    
13. 15  2
    
14. 
    
15. 16  2
    
16. 
    
17. 17  2
    
18. 
    
19. 19  2
    
20. 
    
21. ;
    
22. 
    
23. run;

复制代码

如何用R写个写个函数实现如SAS输出的结果，谢谢！

1. myfunc <- function(..., digits = 2){
   
2.   list <- list(...)
   
3.   if (is.call(list[[1]])) {
   
4.     mydata <- eval(list[[2]])
   
5.     if (!(is.matrix(mydata) || is.data.frame(mydata)))
   
6.       stop("第二个参数须为矩阵或者数据框")
   
7.     listformula <- as.list(list[[1]])
   
8.     x <- mydata[[deparse(listformula[[2]])]]
   
9.     grp <- mydata[[deparse(listformula[[3]])]]
   
10.   } else {
    
11.     x <- c(list[[1]], list[[2]])
    
12.     grp <- rep(c(1, 2), c(length(list[[1]]), length(list[[2]])))
    
13.     mydata <- data.frame(x, grp)
    
14.   }
    
15.   options(digits = digits)
    
16.   
    
17.   satt <- t.test(..., var.equal = F)
    
18.   ttest <- t.test(..., var.equal = T)
    
19.   vartest <- var.test(list[[1]], list[[2]])
    
20.   
    
21.   summ <- by(x, grp, pastecs::stat.desc)
    
22.   lisnames <- names(summ)
    
23.   table1names <- c("g", "方法", "数目", "均值", "标准差", "标准误差",
    
24.                   "最小值", "最大值")
    
25.   table2names <- c("g", "方法", "均值", "%95-down", "%95-up", "标准差")
    
26.   
    
27.   vector1 <- summ[[1]][c(9, 13, 10, 4, 5)]
    
28.   vector2 <- summ[[2]][c(9, 13, 10, 4, 5)]
    
29.   
    
30.   vec1 <- c(summ[[1]][9], summ[[1]][9] - summ[[1]][11],
    
31.             summ[[1]][9] + summ[[1]][11], summ[[1]][13])
    
32.   vec2 <- c(summ[[2]][9], summ[[2]][9] - summ[[2]][11],
    
33.             summ[[2]][9] + summ[[2]][11], summ[[2]][13])
    
34.   
    
35.   diffmean <- summ[[1]][9] - summ[[2]][9]
    
36.   
    
37.   diffsum1 <- c(diffmean, NA, ttest$stderr, NA, NA)
    
38.   diffsat1 <- c(diffmean, NA, satt$stderr, NA, NA)
    
39.   
    
40.   diffsum2 <- c(diffmean, ttest$conf.int[1:2], NA)
    
41.   diffsat2 <- c(diffmean, satt$conf.int[1:2], NA)
    
42.   
    
43.   words <- c("差", "(", lisnames[1], "-", lisnames[2], ")")
    
44.   g <- c(lisnames, paste(words, collapse = ""), paste(words, collapse = ""))
    
45.   method <- c(NA, NA, "汇总", "Satterthwaite")
    
46.   num <- c(length(grp[as.character(grp) == lisnames[1]]),
    
47.            length(grp[as.character(grp) == lisnames[2]]), NA, NA)
    
48.   
    
49.   matrix1 <- rbind(vector1, vector2, diffsum1, diffsat1)
    
50.   table1 <- data.frame(g, method, num, matrix1)
    
51.   colnames(table1) <- table1names
    
52.   rownames(table1) <- NULL
    
53.   
    
54.   matrix2 <- rbind(vec1, vec2, diffsum2, diffsat2)
    
55.   table2 <- data.frame(g, method, matrix2)
    
56.   colnames(table2) <- table2names
    
57.   rownames(table2) <- NULL
    
58.   
    
59.   方法 <- c("汇总", "Satterthwaite")
    
60.   方差 <- c("等于", "不等于")
    
61.   自由度 <- c(ttest[["parameter"]],satt[["parameter"]])
    
62.   t值 <- c(ttest[["statistic"]],satt[["statistic"]])
    
63.   `Pr > |t|` <- c(ttest[["p.value"]],satt[["p.value"]])
    
64.   table3 <- data.frame(方法, 方差, 自由度, t值, `Pr > |t|`)
    
65.   rownames(table3) <- NULL
    
66.   
    
67.   table4 <- data.frame(方法 = "折叠的F",
    
68.                          分子自由度 = vartest$parameter[1],
    
69.                          分母自由度 = vartest$parameter[2],
    
70.                          F值 = vartest$statistic,
    
71.                          `Pr > F` = vartest$p.value)
    
72.   rownames(table4) <- NULL
    
73.   
    
74.   print.ttest <- function(table){
    
75.       print(table)
    
76.       cat("\n")
    
77.   }
    
78.   
    
79.   table <- list(table1, table2, table3, table4)
    
80.   for (i in 1:length(table)) print.ttest(table[[i]])
    
81. }
    
82. 
    
83. #example
    
84. x <- c(12 : 17, 19)
    
85. grp <- rep(2 : 3, c(3, 4))
    
86. testdata <- data.frame(x, grp)
    
87. 
    
88. myfunc(x~grp, testdata)
    
89. # or
    
90. myfunc(x[1:3], x[4:7])
    
91. 
    
92. # digits参数可以控制小数位
    
93. myfunc(x[1:3], x[4:7], digits = 1)
    
94. 
    
95. # alternative可以设置检验方向:c("two.sided", "less", "greater")
    
96. myfunc(x[1:3], x[4:7], alternative = "less")

复制代码

需要安装pastecs包， 两个软件计算方差齐性时的F值略有不同，R是小方差除以大方差,sas是大方差除以小方差,相当于一个是左侧F值，一个右侧F值。有几个描述指标没有找到现成的函数，楼主可以在我的函数里添加这些指标的计算。我就不写了。

方差齐性检验：

1. var.test(x~g, data=)

复制代码

若方差不齐采用satterthwaite法：

1. t.test(x~g, data= , var.equal = F)

复制代码

然后根据把结果整合输出就行了，不是很难。

乐天天12300 发表于 2021-1-17 14:28
方差齐性检验：

若方差不齐采用satterthwaite法：

能否写个函数整合下，结果输出如SAS结果

继续     

继续求助！

继续求助

乐天天12300 发表于 2021-1-16 23:36
需要安装pastecs包， 两个软件计算方差齐性时的F值略有不同，R是小方差除以大方差,sas是大方差除以小方 ...

function(...
函数后面... 什么意思，谢谢！

dxystata 发表于 2021-2-28 16:56
function(...
函数后面... 什么意思，谢谢！

是一种特殊类型的参数。