自定义相同指标内多组数据对比：单因素方差分析函数

注意：无论以哪种形式作为数据输入，符号~右边表示分类变量，符合左边表示因变量，即要研究的指标数据

aov.func<-function(data,alpha=0.05,plot.logic=T,p.adjust.mod="holm"){

+ error<-FALSE

+ data[,1]<-as.factor(data[,1])

+ group.levels<-as.numeric(levels(data[,1]))#分类变量的水平取值

+ for(i in group.levels){

+ if(shapiro.test(data[which(data[,1]==i),2])$p.value<alpha){+ print(paste("ERROR:第",i,"组数据不服从正态分布"))

#验证不同水平下的指标是否服从正态分布，只要有一组不服从正态分布，就带引第i组不服从正态分布

+ error<-TRUE

+ }

+ }

+ if(error){

+ return()

+ }

+ else{

+ print("符合正态性前提")

+ }

+ #检验不同水平下指标的方差齐次性

+ if(bartlett.test(data[,2]~data[,1])$p.value<alpha){

+ print("ERROR:符合方差齐次性前提")

+ }

#验证不同水平下指标数据的方差是否相同，如果p-value大于alpha表示服从，小于则表示不服从

+ #绘制箱线图

+ if(plot.logic){

+ boxplot(data[,2]~data[,1])

+ }

+ #方差分析

+ sol<-aov(data[,2]~data[,1])

+ p.value<-summary(sol)[[1]][,5][1]

+ tab<-matrix(NA,nrow=3,ncol=5)

+ dimnames(tab)<-list(c("分类变量A","随机性误差","总和"),c("自由度","差异平方和","均方","统计量F","P-value"))

+ tab[1:2,1]<-summary(sol)[[1]][,1];tab[3,1]<-sum(tab[1:2,1])#自由度

+ tab[1:2,2]<-summary(sol)[[1]][,2];tab[3,2]<-sum(tab[1:2,2])#差异平方和

+ tab[1:2,3]<-summary(sol)[[1]][,3];#均方

+ tab[1:2,4]<-summary(sol)[[1]][,4];#统计量F

+ tab[1,5]<-p.value;#p-value值

+ print("====方差分析====")

+ print(tab)

+ #多重T检验

+ if(p.value<alpha){

+ print(paste(p.value,"<",alpha,":各水平下的指标有明显差别"))

+ sol.t<-pairwise.t.test(data[,2],data[,1],p.adjust.method=p.adjust.mod)

+ tab<-sol.t[[3]]

+ tab[which((sol.t[[3]]>alpha)==TRUE)]<-"无差别"

+ tab[which((sol.t[[3]]>alpha)==FALSE)]<-"显著差别"

+ print("====多重T检验====")

+ print(tab)

+ }

+ else{

+ print(paste(p.value,">",alpha,":个水平下的指标无明显差别"))

+

+ }

+ return(sol)

+ }

> 自定义函数aov.func的参数说明：

data:可以是matrix或data.frame两种对象，第一列表示不同水平下的取值，第二列是对应的指标数据，例如：kpi<-c(rnorm(20,80,5),rnorm(25,195,5),rnorm(23,256,5)),group<-c(rep(1,20),rep(2,25),rep(3,23))

data<-data.frame(group,kpi),

Alpha表示显著性水平，默认为0.05

Plot.logic如果去T默认，时绘制data的箱线图，如果去F时不进行绘图。

P.adjust.mod设置在左多重T检验时，修正p-value的值。

例：数据为某电商2个年龄段的日客单价数据。

>  price1<-c(79,79,87,79,71,84,82,82,85,82,88,81,71,76,81,75,72,83,76)

>price2<-c(193,192,191,181,191,192,187,191,196,196,196,196,192,194,192,189,196,202,197,206,199,191,194,195,196)

> date<-rep(1:2,c(19,25))

> price<-c(price1,price2)

> data<-data.frame(group=date,x=price)

> data

   group   x

1      1  79

2      1  79

3      1  87

4      1  79

5      1  71

6      1  84

7      1  82

8      1  82

9      1  85

10     1  82

11     1  88

12     1  81

13     1  71

14     1  76

15     1  81

16     1  75

17     1  72

18     1  83

19     1  76

20     2 193

21     2 192

22     2 191

23     2 181

24     2 191

25     2 192

26     2 187

27     2 191

28     2 196

29     2 196

30     2 196

31     2 196

32     2 192

33     2 194

34     2 192

35     2 189

36     2 196

37     2 202

38     2 197

39     2 206

40     2 199

41     2 191

42     2 194

43     2 195

44     2 196

> sol<-aov.func(data)

[1] "符合正态性前提"

[1] "====方差分析===="

           自由度 差异平方和         均方  统计量F      P-value

分类变量A       1 140712.579 140712.57895 5780.327 1.289547e-46

随机性误差     42   1022.421     24.34336       NA           NA

总和           43 141735.000           NA       NA           NA

[1] "1.28954695480505e-46 < 0.05 :各水平下的指标有明显差别"

[1] "====多重T检验===="

  1         

2 "显著差别"

>