数据分析:拯救你的很多时候是基础理论
2023-05-10 14:56:27其实数据挖掘工程师的人力资源一般从两个群体里面而来，其一数学、统计学、应用科学（如生物统计学）相关的学科转行而来，这群人统计学扎实，比较注重分析过程的统计基础和算法，但是代码的速度和简洁经常遭人鄙视，他们一般选择R和matlab之类的分析工具，其二由程序员转行而来，这类人代码质量较高，但是分析过程的统计基础较差，他们一般选择python之类的分析工具。但可笑的是总体而言，这些人都会逐渐忘掉原来所学的统计基础，有时候记住了很多高大上的算法却忘记了基础的统计分析，比如分析哪一种人的收入较高时，或者哪一种基因导致的转氨酶不正常表达时，竟然无法对应到基础的均值比较或方差分析，评价算法的效果时竟然不能自己独立的构建评价指标，有时候即使做了评比，也不具有统计意义，甚至有人滥用“显著”这种具有特殊含义的统计字眼。
4.1.参数检验及非参检验
参数检验最常见的是均值比较，本书牵涉到的统计问题一般是不拿来说一说就无法继续的统计知识，而均值比较就是这一类，它的重要性几乎会纠缠每一个数据挖掘工程师的一生，也会贯穿一个项目周期。“比较”一般是通过参数检验和非参检验进行的，而均值比较就是最常见的参数检验，比如你要比较两个模型的效果，肯定要进行准确率或者mae、mse、rmse、mape等评价指标的均值的比较，单纯的一次建模和测试就进行比较是没有统计意义。
比如我们要比较同时期的人群收入是否存在差异，比如我们要分析三级医院和二级医院的收入差异，这些都要用到均值比较。
均值比较中最常用的是t检验（样本量小于30用z检验），t检验要求数据符合正态分布的假设，所以在进行t检验之前首先要进行正态分布的检验，正态分布的检验分为两类，一种是公式法，如Shapiro-Wilk、Anderson-Darling、Kolmogorov-Smirnov等检验方法，另外一种是图形法，如qq图、频度直方图等等，作为一个经常要求分析自动化进行的业内人士，我更倾向于第一种，符合正态分布才能进行t检验，但是这个两步法并非没有问题，详情参看《Is normality testing ‘essentially useless’?》。
首先进行一个单样本的均值检验，比如我们想知道北京的二级医院平均收入是否高于1亿元，这就用到了均值检验，这类问题仅包含一个样本和具体值的比较，所以称为单样本t检验。首先要检验北京的医院收入是否符合正态分布，然后再进行t检验。
qq图检验正态分布
hospital <- read.csv("H:/zimeiti/窥视数据背后的逻辑：基于R与python/bookwriting/第二章数据探索招招都是利器/data/hospital.csv",
header = T, sep = ",", stringsAsFactors = F)
bjhospital <- hospital[hospital$省份 == "北京" & hospital$医院级别 == "二级",]
summary(bjhospital$总收入)
# Min. 1st Qu. Median Mean 3rd Qu. Max.
# 0 27250 84790 125500 181200 548200
qqnorm(bjhospital$总收入)
qqline(bjhospital$总收入, col = 2)
# 标准的正态分布
temp <- rnorm(1000, mean = 2)
qqnorm(temp)
qqline(temp, col = 2)图形检验正态分布通常在ppt展示时比较常用，方法有频率图、密度图、qq图，它们有共同的优点比较直观，1行读取医院的数据；2行筛选出北京市的二级医院，“&”在R里表示且的关系，3行使用summary看一下北京二级医院总收入的统计摘要，可以看到北京医院的平均收入为1.2亿，那我们现在是否可以说北京二级医院总收入平均显著大于1亿元呢？显然这个时候做出这样的表述有点不负责任，且显得缺乏专业素养，因为还没进行统计意义上的比较；4行qqnorm绘制qq图，qq是quantile的缩写，用于检验样本的理论百分位数分布于实际百分位数是否相符，如果相符，纵轴是数据的样本分位数，横轴是正态分布的理论分位数，沿45度的参考线；qqline函数为qq图添加标准的拟合线，如果数据是正态分布，则散点应该落在直线上，即实际分位数与理论分位数之间偏差很小，反之，如果偏离较大，则说明数据并非来自正态总体，但是我们看到北京二级医院的总收入并没有与理论上完全相符，显然有些扭曲，至少不是标准的正态分布，这时图形法的缺点就出现了，即依据图形法判断实在太过主观，一旦和标准不太一样时，只能根据个人经验判断；7-9行我们我们绘制了一张标准的qq图，首先我们使用rnorm产生了一个均值为2的标准分布，8行绘制qq图，9行添加拟合线，我们看到标准的正态分布基本上沿线分布，看了这幅图等于给自己一个标准的印象。
正态分布检验
shapiro.test(bjhospital$总收入)
# Shapiro-Wilk normality test
#
# data: bjhospital$总收入
# W = 0.84748, p-value = 2.334e-08
ks.test(bjhospital$总收入, "pnorm")
# One-sample Kolmogorov-Smirnov test
#
# data: bjhospital$总收入
# D = 0.98925, p-value < 2.2e-16
# alternative hypothesis: two-sided常见的检验正态分布的方法之前说过，但是我个人更加倾向于选择Kolmogorov-Smirnov方法，Shapiro-Wilk方法对样本大小的依赖性很大，一般样本量小于50，这也是shapiro.test函数将样本量限制在5000以内的原因，如果样本比较大，Shapiro-Wilk很容易拒绝原假设得出不符合正态分布的结论，比如1行使用shapiro.test检查总收入的正态性，看到p-value远远小于0.05甚至远远低于0.01，所以拒绝原假设，北京二级医院的总收入不符合正态分布，如果p-value远远大于0.05表示不能拒绝原假设，也就是说不能够否定数据属于正态分布，很遗憾北京二级医院的收入不符合正态分布，因为ks.test函数的检验结果也支持这种结论。
但是也不要太担心，其实有时候在实战操作中如果样本量很大，就直接认为符合正态分布了，另外，即使原假设被拒绝了，也可以根据数据背景讨论，数据可能受某些异常值的影响，偏离了正态分布，可以根据上面讲到分析异常值的方法将异常值剔除之后再检验，实在通不过还可以进行非参检验，这里我们先假设数据符合正态分布，继续往下一步走。
单样本t检验
t.test(x = bjhospital$总收入, alternative = "greater", mu = 100000)
# One Sample t-test
#
# data: bjhospital$总收入
# t = 1.9485, df = 92, p-value = 0.0272
# alternative hypothesis: true mean is greater than 1e+05
# 95 percent confidence interval:
# 103762 Inf
# sample estimates:
# mean of x
# 125546.71行使用t.test函数进行均值比较，比较北京二级医院收入的均值是否大于1亿，第一个参数用于指定样本数据，alternative指定检验类型，这里指定单尾检验，我们看到t值为1.9，p-value小于0.05，我们可以认为拒绝原假设，认为北京二级医院平均收入显著大于1亿。
除了单样本的均值检验以外，我们还会遇到双样本的均值检验，在进行t检验之前一般认为两个样本来自于方差相同的正态总体，比如我们检验北京医院的平均收入和河南医院的平均收入。
双样本t检验
bjhospital <- hospital[hospital$省份 == "北京", "总收入"]
hnhospital <- hospital[hospital$省份 == "河南", "总收入"]
shapiro.test(bjhospital)
# Shapiro-Wilk normality test
#
# data: bjhospital
# W = 0.40661, p-value < 2.2e-16
shapiro.test(hnhospital)
# Shapiro-Wilk normality test
#
# data: hnhospital
# W = 0.29183, p-value < 2.2e-16
result <- t.test(x = bjhospital, y = hnhospital, var.equal = FALSE)
# Welch Two Sample t-test
#
# data: bjhospital and hnhospital
# t = 6.0533, df = 579.38, p-value = 2.55e-09
# alternative hypothesis: true difference in means is not equal to 0
# 95 percent confidence interval:
# 39078.56 76617.38
# sample estimates:
# mean of x mean of y
# 79660.77 21812.80
result
result$statistic
result$p.value
result$conf.int1行提取北京医院的的数据；2行提取河南医院的数据，3-4行检验样本是否符合正态分布，结果两个省市的医院收入都不符合正态分布，按道理下面就不能轻易进行t检验，但是我们这里为了演示，继续进行；5行进行t检验，同时设置了方差不相等，如果你的数据方差相等就将此项设置为T，另外需要注意的是t.test还有一个非常重要的参数，paired = FALSE，用于指定是否进行配对数据检验，如果设为T，表示数据是一一配对的，比如对同一批病人一个1月转氨酶与2月转氨酶是否存在差异，就可以使用配对数据检验，因为检验的是同一个样本的先后差异；6行查看结果，看到北京医院的平均收入为7966万，河南的医院平均收入为2181万，他们之间差值的95%的置信区间为[39078.56,76617.38],因为p-value远远小于0.01，所以拒绝原假设，认为两个省份的医院收入存在极显著差异；有时候我们需要提取结果的一些统计量存储备用，比如7-9行粉笔提取了t值、pvalue和差异置信区间。
除了均值检验以外，我们可能还会进行非参检验，非参检验一般是在数据不符合参数检验要求时进行，比如上一步我们的数据不符合正态分布的要求，我们就可以进行非参检验，这里我们重点演示一下秩和检验和卡方检验。
秩和检验适用于上一步中检验不符合正态分布的情况，比如我们比较北京二级医院平均收入是否大于1亿，但样本不符合正态分布，我们就可以使用秩和检验，秩和检验是进行中位数之间的比较。
单样本秩和检验
bjhospital <- hospital[hospital$省份 == "北京" & hospital$医院级别 == "二级",]
summary(bjhospital$总收入)
# Min. 1st Qu. Median Mean 3rd Qu. Max.
# 0 27250 84790 125500 181200 548200
wilcox.test(bjhospital$总收入, alternative = "greater", mu = 100000)
# Wilcoxon signed rank test with continuity correction
#
# data: bjhospital$总收入
# V = 2297, p-value = 0.3353
# alternative hypothesis: true location is greater than 1e+051行提取北京市二级医院的收入；2行查看统计摘要，可以看到北京市二级医院的收入中位数为8500万左右，是不是存在大于1亿的可能呢？3行我们使用wilcox.test函数进行了wilcox秩和检验，发现结果p-value等于0.33，远大于0.05，所以不认为北京二级医院收入的中位数大于1亿元。
秩和检验
bjhospital <- hospital[hospital$省份 == "北京", "总收入"]
hnhospital <- hospital[hospital$省份 == "河南", "总收入"]
wilcox.test(x = bjhospital, y = hnhospital)
# Wilcoxon rank sum test with continuity correction
#
# data: bjhospital and hnhospital
# W = 347210, p-value = 0.0005514
# alternative hypothesis: true location shift is not equal to 0
wilcox.test(x = bjhospital, y = hnhospital, alternative = "greater")
# Wilcoxon rank sum test with continuity correction
#
# data: bjhospital and hnhospital
# W = 347210, p-value = 0.0002757
# alternative hypothesis: true location shift is greater than 01-2提取北京和河南的医院收入，3行wilcox.test函数检验两省市的医院收入中位数是否存在差异，p-value远远小于0.01，更不用说95%的置信度，可以确认两省市的医院收入中位数存在差异；4行wilcox.test函数检验北京市医院收入的中位数是否大于河南省的，结果p-value远远小于0.01，接受备择假设认为北京医院收入的中位数极显著大于河南省医院收入中位数。
另外我们经常需要检验两个变量是否相关，比如医院收入是否和医院的床位数相关，这是就用到了cor.test函数，用于比较相关性。
相关性检验
cor.test(x = hospital$总收入, y = hospital$床位数, method = "pearson")
# Pearson's product-moment correlation
#
# data: hospital$总收入 and hospital$床位数
# t = 190.48, df = 19169, p-value < 2.2e-16
# alternative hypothesis: true correlation is not equal to 0
# 95 percent confidence interval:
# 0.8039421 0.8137304
# sample estimates:
# cor
# 0.8088923cor.test函数用于计算两个变量的相关性，结果发现总收入和床位数存在非常明显的正相关，相关系数为0.809，p-value远远小于0.01，接受备择假设，认为两个变量之间存在极显著的相关性，cor.test函数参数method用于指定计算相关系数的方法，如果是连续性变量用pearson，如果是等级相关的变量请用spearman，即秩相关，比如比较两个排名结果是否存在差异等等，就可以转化成两个排名的spearman相关系数的检验，排名章节将详细介绍，与之相比，我们更需要尽快了解卡方检验。比如我们要检验样本是否符合某种概率分布，就要用到卡方检验了。
例如我们分析大麦杂交后形状比例是否符合无芒：长芒：短芒=9:3:4，实际观测值435：85：146；则可以进行如下检验。
卡方检验
prob <- c(9/16,3/16,4/16)
x <- c(435, 85, 146)
chisq.test(x, p = prob)
# Chi-squared test for given probabilities
#
# data: x
# X-squared = 24.987, df = 2, p-value = 3.751e-061生成理论概率，2创建对象x，按照概率类别的顺序输入数值；3行chisq.test函数进行卡方检验，结果p-value远远小于0.01，拒绝原假设，认为真实值和概率分布存在非常明显的差异。
另外卡方检验通常用于检验两个现象是否是独立的，比如吸烟和患肺癌是否是独立现象，如下例：
卡方独立性检验
x <- matrix(c(90,10,12,31),nrow=2)
colnames(x) <- c("患癌", "未患癌")
rownames(x) <- c("吸烟", "不吸烟")
chisq.test(x)
# Fisher's Exact Test for Count Data
#
# data: x
# p-value = 2.31e-13
# alternative hypothesis: true odds ratio is not equal to 1
# 95 percent confidence interval:
# 8.386529 66.055954
# sample estimates:
# odds ratio
# 22.44285
fisher.test(x)
# Fisher's Exact Test for Count Data
#
# data: x
# p-value = 2.31e-13
# alternative hypothesis: true odds ratio is not equal to 1
# 95 percent confidence interval:
# 8.386529 66.055954
# sample estimates:
# odds ratio
# 22.442851行创建样本统计矩阵；2行修改矩阵列名；3行修改矩阵行名称；4行chisq.test函数进行卡方独立性检验，结果p-value远远小于0.01，认为二者之间存在非常强烈的联系；有时chisq检验的敏感度较低，可以使用费雪检验（fisher.test），结果同样表明二者之间存在非常明显的相关性。
基本上到这里我们将统计学的基础知识参数检验和非参检验过了一遍，本书并非完全的统计教参，疏漏难免，但是可以肯定的是数据挖掘之道就是要遵循基本的统计学，脱离了统计逻辑的报告和过程都是离经叛道之举，陷自己于不专业的境地，比如没有进行合格的统计比较就乱用“显著”和“极显著”这样的字眼等等。
链接: https://pan.baidu.com/s/1o7Cex6a 密码: 4stc