如何解决R:用于正态性检验的 Kolmogorov-Smirnov 检验可视化经验 CDF 与理论 CDF 的比较
我想使用 Kolmogorov-Smirnov 检验 (ks-test) 检查正态性并将其可视化以支持我的演示。
我使用了这个默认代码
sample1 <- rnorm(10000,10,5)
sample2 <- rnorm(10000,1,5)
group <- c(rep("sample1",length(sample1)),rep("sample2",length(sample2)))
dat <- data.frame(KSD = c(sample1,sample2),group = group)
# create ECDF of data
cdf1 <- ecdf(sample1)
cdf2 <- ecdf(sample2)
# find min and max statistics to draw line between points of greatest distance
minMax <- seq(min(sample1,max(sample1,length.out=length(sample1))
x0 <- minMax[which( abs(cdf1(minMax) - cdf2(minMax)) == max(abs(cdf1(minMax) - cdf2(minMax))) )]
y0 <- cdf1(x0)
y1 <- cdf2(x0)
ggplot(dat,aes(x = KSD,group = group,color = group))+
stat_ecdf(size=1) +
theme_bw(base_size = 28) +
theme(legend.position ="top") +
xlab("Sample") +
ylab("ECDF") +
#geom_line(size=1) +
geom_segment(aes(x = x0[1],y = y0[1],xend = x0[1],yend = y1[1]),linetype = "dashed",color = "red") +
geom_point(aes(x = x0[1],y= y0[1]),color="red",size=8) +
geom_point(aes(x = x0[1],y= y1[1]),size=8) +
ggtitle("K-S Test: Sample 1 / Sample 2") +
theme(legend.title=element_blank())
现在,我想用我自己的数据替换它
mydata <- c(33,28,23,27,22,25,32,26,14,37,16,30,34,18,31,20,24,12,40,17,29,36,46,33,21,35,13,19,39,11,15,42,38,50,44,41,27)
我可以将 mydata
用作经验 CDF,但如何将 sample2
更改为理论 CDF?
目标:
- 2 cdf 图:经验 CDF(
mydata
替代sample1
)和理论 CDF(替代sample2
分布 - 出现最大差异的垂直线 +
geom_pint
用于加法
提前致谢。
解决方法
R 中的理论累积概率函数由分布函数的 p*
版本给出。如果你想测试正态性,函数应该是 pnorm
。在下面,您将找到一种绘制此图的方法。黑线表示最大差异点。
library(ggplot2)
#> Warning: package 'ggplot2' was built under R version 4.0.5
mydata <- c(33,28,23,10,27,22,25,32,26,14,37,16,30,34,18,31,20,24,12,40,17,29,36,46,33,21,35,13,19,39,11,15,42,38,50,44,41,27)
# Find x position of maximum difference.
# You could precalculate `mean(mydata)` and `sd(mydata)` if you wish,# we use it later in plotting too.
maxdiff <- mydata[which.max(
abs(ecdf(mydata)(mydata) - pnorm(mydata,mean(mydata),sd = sd(mydata)))
)]
df <- data.frame(x = mydata)
ggplot(df,aes(x)) +
stat_ecdf(aes(colour = "Empirical")) +
stat_function(fun = pnorm,args = list(mean = mean(mydata),sd = sd(mydata)),aes(colour = "Theoretical")) +
annotate(
"segment",x = maxdiff,xend = maxdiff,y = ecdf(mydata)(maxdiff),yend = pnorm(maxdiff,sd(mydata)),size = 1
)
由 reprex package (v1.0.0) 于 2021 年 4 月 25 日创建
版权声明:本文内容由互联网用户自发贡献,该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 dio@foxmail.com 举报,一经查实,本站将立刻删除。