原文：annas-archive.org/md5/60a7ce9fba25ab267b2031ea2a27f285

译者：飞龙

协议：CC BY-NC-SA 4.0

前言

统计学是一门应用分析方法的学科，用于在学术和行业环境中使用数据来回答问题和解决问题。许多方法已经存在了几个世纪，而其他方法则更为新颖。统计分析和结果对于向技术和非技术受众展示来说相当直接。此外，使用统计分析产生结果并不一定需要大量数据或计算资源，而且可以相当快速地完成，尤其是在使用 Python 等易于处理和实现的编程语言时。

随着计算能力的提升和可访问性的增加，近年来人工智能（AI）和高级机器学习（ML）工具变得更加突出和受欢迎。在进行统计分析作为使用 AI 和 ML 开发更大规模项目的前置条件时，这可以使从业者在使用更大计算资源和项目架构开发之前评估可行性和实用性。

本书提供了一系列常用的工具，这些工具可以用于测试假设并为分析师和数据科学家提供基本的预测能力。在探索不同测试及其适用条件之前，读者将了解理解本书中统计工具所需的基本概念和术语。此外，读者还将获得评估测试性能的知识。在整个过程中，将通过 Python 编程语言提供示例，帮助读者开始使用所提供的工具理解他们的数据，这些工具将适用于数据分析行业面临的一些最常见问题。我们将探讨以下主题：

• 统计学简介

• 回归模型

• 分类模型

• 时间序列模型

• 生存分析

理解这些部分提供的工具将为读者提供一个坚实的基础，从而更容易在统计学领域实现进一步的独立增长。

本书面向的对象

大多数行业的专业人士都可以从本书中的工具中受益。提供的工具主要用于更高层次的推理分析，但根据从业者希望应用的行业，也可以应用于更深的层次。本书的目标受众包括：

• 拥有限统计或编程知识但希望学习如何使用数据来测试他们在商业领域持有的假设的行业专业人士

• 希望扩展他们的统计知识并找到用于执行各种数据导向任务的工具及其实现的数据分析师和科学家

本书自下而上的方法旨在为广泛的受众提供进入知识库的途径，因此不应使新手级从业者感到气馁，也不应排除高级从业者从所提供材料中获益。

本书涵盖的内容

第一章, 抽样与泛化，描述了抽样和泛化的概念。抽样讨论涵盖了从总体中抽取数据的几种常见方法，并讨论了泛化的影响。本章还讨论了如何设置本书所需的软件。

第二章, 数据分布，详细介绍了数据类型、描述数据的常用分布和统计量。本章还涵盖了用于改变分布的常见转换。

第三章, 假设检验，介绍了统计检验的概念，作为回答感兴趣问题的方法。本章涵盖了进行检验的步骤、测试中遇到的错误类型，以及如何使用 Z 检验选择功效。

第四章, 参数检验，进一步讨论了统计检验，提供了常见参数统计检验的详细描述、参数检验的假设以及如何评估参数检验的有效性。本章还介绍了多重检验的概念，并提供了多重检验校正的详细信息。

第五章, 非参数检验，讨论了当参数检验的假设被违反时如何进行统计检验，即使用无假设的测试类别，称为非参数检验。

第六章, 简单线性回归，通过简单线性回归模型介绍了统计模型的概念。本章首先讨论简单线性回归的理论基础，然后讨论如何解释模型的结果以及评估模型的有效性。

第七章, 多元线性回归，在前一章的基础上，将简单线性回归模型扩展到更多维度。本章还讨论了使用多个解释变量建模时出现的问题，包括多重共线性、特征选择和降维。

第八章, 离散模型，介绍了分类的概念，并开发了一个将变量分类到分类响应变量的离散水平中的模型。本章首先开发二元分类模型，然后将模型扩展到多元分类。最后，介绍了泊松模型和负二项式模型。

第九章，判别分析，讨论了几个额外的分类模型，包括线性判别分析和二次判别分析。本章还介绍了贝叶斯定理。

第十章，时间序列简介，介绍了时间序列数据，讨论了时间序列的自相关概念和时间序列的统计度量。本章还介绍了白噪声模型和稳定性。

第十一章，ARIMA 模型，讨论了单变量模型。本章首先讨论了平稳时间序列的模型，然后将讨论扩展到非平稳时间序列。最后，本章提供了模型评估的详细讨论。

第十二章，多元时间序列，在前两章的基础上介绍了多元时间序列的概念，并将 ARIMA 模型扩展到多个解释变量。本章还讨论了时间序列的交叉相关性。

第十三章，生存分析，介绍了生存数据，也称为事件时间数据。本章讨论了截尾的概念和截尾对生存数据的影响。最后，本章讨论了生存函数、风险率和风险比。

第十四章，生存模型，在前一章的基础上，概述了几个生存数据模型，包括 Kaplan-Meier 模型、指数模型和 Cox 比例风险模型。

为了充分利用本书

您需要访问下载和安装 Python 编程语言中实现的开放源代码包，这些包可通过 PyPi.org 或 Anaconda Python 发行版访问。虽然统计学背景有帮助，但不是必需的，本书假设您有扎实的代数基础知识。本书的每个单元都是独立的，但每个单元内的章节是相互关联的。因此，我们建议您从每个单元的第一章开始，以了解 内容 。

本书涵盖的软件/硬件	操作系统要求
Python 版本 ≥ 3.8	Windows, macOS, 或 Linux
Statsmodels 0.13.2
SciPy 1.8.1
lifelines 0.27.4
scikit-learn 1.1.1
pmdarima 2.02
Sktime 0.15.0
Pandas 1.4.3
Matplotlib 3.5.2
Numpy 1.23.0

如果您正在使用本书的数字版，我们建议您亲自输入代码或从本书的 GitHub 仓库（下一节中提供链接）获取代码。这样做将有助于避免与代码复制和粘贴相关的任何潜在错误。

下载示例代码文件

您可以从 GitHub 下载此书的示例代码文件github.com/PacktPublishing/Building-Statistical-Models-in-Python。如果代码有更新，它将在 GitHub 仓库中更新。

我们还有其他来自我们丰富的图书和视频目录的代码包可供选择，可在github.com/PacktPublishing/找到。查看它们！

使用的约定

本书使用了多种文本约定。

文本中的代码：表示文本中的代码单词、数据库表名、文件夹名、文件名、文件扩展名、路径名、虚拟 URL、用户输入和 Twitter 昵称。以下是一个示例：“将下载的WebStorm-10*.dmg磁盘映像文件作为系统中的另一个磁盘挂载。”

代码块设置如下：

A = [3,5,4]
B = [43,41,56,78,54]
permutation_testing(A,B,n_iter=10000)

任何命令行输入或输出都按照以下方式编写：

pip install SomePackage

粗体：表示新术语、重要单词或你在屏幕上看到的单词。例如，菜单或对话框中的单词以粗体显示。以下是一个示例：“从管理面板中选择系统信息。”

小贴士或重要注意事项

看起来是这样的。

联系我们

我们始终欢迎读者的反馈。

一般反馈：如果您对此书的任何方面有疑问，请通过 customercare@packtpub.com 给我们发邮件，并在邮件主题中提及书名。

勘误：尽管我们已经尽一切努力确保内容的准确性，但错误仍然可能发生。如果您在此书中发现错误，我们将不胜感激，如果您能向我们报告，我们将不胜感激。请访问www.packtpub.com/support/errata并填写表格。

盗版：如果您在互联网上以任何形式遇到我们作品的非法副本，如果您能提供位置地址或网站名称，我们将不胜感激。请通过 copyright@packtpub.com 与我们联系，并提供材料的链接。

如果你有兴趣成为作者：如果你在某个领域有专业知识，并且你感兴趣的是撰写或为书籍做出贡献，请访问authors.packtpub.com。

分享你的想法

读完《使用 Python 构建统计模型》后，我们很乐意听听您的想法！请点击此处直接进入此书的亚马逊评论页面并分享您的反馈。

您的评论对我们和科技社区都很重要，并将帮助我们确保我们提供高质量的内容。

下载此书的免费 PDF 副本

感谢您购买此书！

你喜欢在路上阅读，但无法携带你的印刷书籍到处走？你的电子书购买是否与你的选择设备不兼容？

别担心，现在每购买一本 Packt 图书，你都可以免费获得该书的 DRM 免费 PDF 版本。

在任何地方、任何设备上阅读。直接从您喜欢的技术书籍中搜索、复制和粘贴代码到您的应用程序中。

优惠不止于此，您还可以获得独家折扣、时事通讯和每日免费内容的每日邮箱访问权限

按照以下简单步骤获取好处：

1. 扫描下面的二维码或访问以下链接

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_QR_Free_PDF.jpg

packt.link/free-ebook/978-1-80461-428-0

1. 提交您的购买证明

2. 就这样！我们将直接将您的免费 PDF 和其他优惠发送到您的电子邮件

第一部分：统计学简介

本部分将涵盖统计建模的基础统计概念。

它包括以下章节：

• 第一章，抽样与推广

• 第二章，数据分布

• 第三章，假设检验

• 第四章，参数检验

• 第五章，非参数检验

第一章：采样和泛化

在本章中，我们将描述种群的概念以及从种群中进行采样的方法，包括一些常见的采样策略。关于采样的讨论将引出描述泛化的部分。泛化将讨论其与使用样本对其各自种群做出结论的关系。在统计推断建模时，确保样本可以推广到种群是必要的。我们将通过本章的主题深入概述这一桥梁。

我们将涵盖以下主要主题：

• 软件和环境设置

• 总体与样本

• 从样本中进行总体推断

• 采样策略 – 随机、系统性和分层

软件和环境设置

Python 是数据科学和机器学习中最受欢迎的编程语言之一，这得益于推动这些库发展的庞大开源社区。Python 的易用性和灵活性使其成为数据科学领域的首选语言，在数据科学领域，实验和迭代是开发周期的关键特征。尽管有新的语言正在开发中用于数据科学应用，例如 Julia，但由于其广泛的开源项目，支持从统计建模到深度学习等应用，Python 目前仍然是数据科学的关键语言。我们选择在本书中使用 Python，因为它在数据科学中的地位重要，并且在就业市场上需求量大。

Python 可用于所有主要操作系统：Microsoft Windows、macOS 和 Linux。此外，安装程序和文档可在官方网站找到：www.python.org/。

本书是为 Python 3.8 版本（或更高版本）编写的。建议您使用可用的最新 Python 版本。本书中的代码不太可能与 Python 2.7 兼容，并且大多数活跃的库已经从 2020 年官方支持结束以来开始停止对 Python 2.7 的支持。

本书使用的库可以使用 Python 包管理器 pip 安装，它是当代 Python 标准库的一部分。有关 pip 的更多信息，请参阅此处：docs.python.org/3/installing/index.xhtml。安装 pip 后，可以在命令行中使用 pip 安装包。以下是一些基本用法：

使用最新版本安装新包：

pip install SomePackage

使用特定版本安装包，例如本例中的版本 2.1：

pip install SomePackage==2.1

已安装的包可以使用 --upgrade 标志进行升级：

pip install SomePackage –upgrade

通常情况下，建议在项目之间使用 Python 虚拟环境，并将项目依赖项与系统目录分开。Python 提供了一个虚拟环境工具 venv，它像 pip 一样是 Python 当代版本标准库的一部分。虚拟环境允许您创建独立的 Python 二进制文件，其中每个 Python 二进制文件都有其自己的安装依赖项集。使用虚拟环境可以防止在处理多个 Python 项目时出现包版本问题和冲突。有关设置和使用虚拟环境的详细信息，请参阅：docs.python.org/3/library/venv.xhtml。

虽然我们推荐使用 Python 和 Python 的虚拟环境进行环境设置，但一个高度推荐的替代方案是 Anaconda。Anaconda 是由 Anaconda Inc.（之前为 Continuum Analytics）提供的免费（企业级）分析型 Python 发行版。Anaconda 发行版包含了许多核心数据科学包、常见的 IDE（如 Jupyter 和 Visual Studio Code）以及用于管理环境的图形用户界面。Anaconda 可以通过在 Anaconda 网站上找到的安装程序进行安装：www.anaconda.com/products/distribution。

Anaconda 自带其自己的包管理器 conda，它可以像 pip 一样用于安装新包。

使用最新版本安装新包：

conda install SomePackage

升级已安装的包：

conda upgrade SomePackage

在本书的整个过程中，我们将使用 Python 数据科学生态系统中的一些核心库，例如 NumPy 用于数组操作，pandas 用于高级数据操作，以及 matplotlib 用于数据可视化。本书使用的包版本包含在以下列表中。请确保您环境中安装的版本等于或高于列表中列出的版本。这将有助于确保代码示例能够正确运行：

• statsmodels 0.13.2

• Matplotlib 3.5.2

• NumPy 1.23.0

• SciPy 1.8.1

• scikit-learn 1.1.1

• pandas 1.4.3

本书代码中使用的包在此处以 图 1*.1* 的形式展示。可以使用 __version__ 方法在代码中打印包版本。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_01_001.jpg

图 1.1 – 本书使用的包版本

在为本书搭建好技术环境之后，让我们进入统计学的讨论。在接下来的章节中，我们将讨论总体和抽样的概念。我们将通过代码实现来演示抽样策略。

总体与样本

通常，统计模型的目标是通过对该群体进行推断来回答关于该群体的问题。我们进行推断的群体可能是生产工厂中的机器、选举中投票的人，或者不同地块上的植物。整个群体，每个单独的项目或实体，被称为总体。在大多数情况下，感兴趣的总体非常大，以至于收集总体中每个实体的数据既不实际也不可能。例如，使用投票的例子，可能无法对选举中投票的每个人进行调查。即使有可能接触到感兴趣的选举的所有选民，许多选民可能不会同意调查，这将阻止对整个总体的收集。另一个考虑因素是调查如此大群体的费用。这些因素使得在我们的投票调查示例中收集总体统计数据实际上是不可能的。在我们可能想要评估总体属性的情况下，存在许多这类禁止性因素。幸运的是，我们不需要收集感兴趣总体中的所有数据。可以使用总体的一部分来进行关于总体的推断。这个总体的子集被称为样本。这是统计模型的主要思想。将使用样本创建模型，并对总体进行推断。

为了使用样本对感兴趣的总体进行有效的推断，样本必须代表感兴趣的总体，这意味着样本应该包含在总体中发现的变异。例如，如果我们对田野中的植物感兴趣，那么来自田野一个角落的样本可能不足以对更大的总体进行推断。整个田野中植物特征可能会有所不同。我们可以考虑各种可能存在变异的原因。对于这个例子，我们将考虑来自图 1.2的一些例子。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_01_002.jpg

图 1.2 – 植物田地

该图显示样本 A靠近森林。这个样本区域可能受到森林存在的影响；例如，该样本中的一些植物可能比其他样本中的植物接收到的阳光更少。样本 B显示位于主要灌溉线之间。可以想象这个样本平均接收到的水可能比其他两个样本多，这可能会影响这个样本中的植物。最后的样本 C靠近道路。这个样本可能看到在样本 A或样本 B中看不到的其他影响。

如果样本只从这些部分中的一个抽取，那么从这些样本中得出的推断将会是有偏的，并且不会提供关于人群的有效参考。因此，需要从整个领域抽取样本，以创建一个更有可能代表植物人群的样本。在从人群抽取样本时，确保抽样方法对可能的问题具有鲁棒性，例如前例中灌溉和遮荫的影响，是至关重要的。每次从人群抽取样本时，都重要的是要识别和减轻可能的偏差影响，因为数据中的偏差会影响你的模型并歪曲你的结论。

在下一节中，将讨论从数据集中进行抽样的各种方法。另一个考虑因素是样本大小。样本大小影响我们可以使用的统计工具类型，关于样本的分布假设，以及推断和预测的置信度。样本大小的影响将在第二章**，数据分布和*第三章**，假设检验中深入探讨。

从样本中进行人群推断

当使用统计模型从该人群的样本子集中对人群进行推断性结论时，研究设计必须考虑到其变量与人群中的变量具有相似的不确定性程度。这就是本章前面提到的变异。为了适当地对人群进行推断性结论，任何统计模型都必须围绕一个随机机制来构建。围绕这些随机机制构建的研究被称为随机实验，并提供了对相关性和因果关系的理解。

随机实验

随机实验有两个主要特征：

• 随机抽样，俗称随机选择

• 治疗的随机分配，这是研究的本质

随机抽样

随机抽样（也称为随机选择）旨在创建一个代表总体人口的样本，以便统计模型能够足够好地推广总体，从而分配因果关系的结果。为了使随机抽样成功，感兴趣的总体必须定义良好。从总体中抽取的所有样本都必须有被选中的机会。在考虑民意调查选民的情况时，所有选民都必须愿意接受调查。一旦所有选民都进入彩票，就可以使用随机抽样来对选民进行子集建模。仅从愿意接受调查的选民中进行抽样会引入统计模型中的抽样偏差，这可能导致结果偏斜。只有一些选民愿意参与的情景中的抽样方法被称为自我选择。从自我选择的样本——或任何非随机样本——中获得并建模的信息不能用于推断。

治疗的随机分配

治疗的随机分配涉及两个动机：

• 第一个动机是了解特定的输入变量及其对响应的影响——例如，了解将治疗 A 分配给特定个体是否可能比安慰剂产生更理想的结果。

• 第二个动机是消除外部变量对研究结果的影响。这些被称为混杂变量（或混杂因素）的外部变量很重要，因为它们往往难以控制。它们可能具有不可预测的值，甚至可能对研究者来说是未知的。包括混杂变量的后果是，研究的结果可能无法复制，这可能是代价高昂的。虽然混杂变量可以影响结果，但它们也可以影响输入变量，以及这些变量之间的关系。

回顾前一小节中的例子，即总体与样本，考虑一位农民决定开始在他的作物上使用杀虫剂，并想要测试两种不同的品牌。这位农民知道土地上有三个不同的区域；地块 A、地块 B 和地块 C。为了确定杀虫剂的成功率并防止作物受损，农民从每个地块随机选择 60 株植物进行测试（这被称为分层随机抽样，即在每个地块上进行分层随机抽样）。这个选择代表了植物总体。从这个选择中，农民对他的植物进行标记（标记不需要是随机的）。对于每个地块，农民将标签随机放入袋子中，以随机化它们，然后开始选择 30 株植物。前 30 株植物接受一种处理，另外 30 株接受另一种处理。这是随机分配的处理。假设三个单独的地块代表了一组不同的混杂变量对作物产量的影响，农民将拥有足够的信息来对每个杀虫剂品牌的作物产量进行推断。

观察性研究

常见的另一种统计研究类型是观察性研究，其中研究人员通过观察已有的数据来寻求了解。观察性研究有助于理解输入变量及其与目标和彼此之间的关系，但不能像随机实验那样提供因果关系理解。观察性研究可能具有随机实验的两个组成部分之一——要么是随机抽样，要么是随机分配的处理——但没有这两个组成部分，将不会直接产生推断。有许多原因可能导致进行观察性研究而不是随机实验，例如以下情况：

• 随机实验成本过高

• 实验的伦理限制（例如，确定怀孕期间吸烟导致出生缺陷发生率的实验）

• 使用先前随机实验的数据，从而消除了进行另一项实验的需要

从观察性研究中推导出一些因果关系的 一种方法是进行随机抽样和重复分析。重复随机抽样和分析可以帮助最小化混杂变量随时间的影响。这个概念在大数据和机器学习的有用性中扮演着巨大的角色，这两者在本世纪许多行业中都获得了很大的重视。虽然几乎任何可以用于观察性分析的工具也可以用于随机实验，但本书主要关注的是观察性分析的工具，因为在大多数行业中这更为常见。

可以说，统计学是一门帮助在存在可量化不确定性时做出最佳决策的科学。所有统计检验都包含一个零假设和一个备择假设。也就是说，一个假设数据之间没有统计学上显著的差异（零假设）或者数据之间存在统计学上显著的差异（备择假设）。统计学上显著差异这一术语意味着存在一个基准——或阈值——超过这个基准，一个度量就会发生并表明其显著性。这个基准被称为临界值。

应用到这个临界值上的度量被称为检验统计量。临界值是基于数据中的行为（如平均值和变异）量化的静态值，并基于假设。如果有两条可能的路径可以拒绝零假设——例如，我们相信某些输出要么小于平均值，要么大于平均值——将有两个临界值（这种检验称为双尾假设检验），但如果只有一条反对零假设的论据，则只有一个临界值（这称为单尾假设检验）。无论临界值的数量有多少，在给定的假设检验中，每个组内都只有一个检验统计量度量。如果检验统计量超过临界值，就有统计学上显著的依据来支持拒绝零假设，并得出数据之间存在统计学上显著差异的结论。

有必要理解假设检验可以检验以下内容：

• 一个变量对另一个变量（例如，在 t 检验中）

• 多个变量对一个变量（例如，线性回归）

• 多个变量对多个变量（例如，多变量方差分析）

在下面的图中，我们可以直观地看到检验统计量与双尾假设检验中的临界值之间的关系。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_01_003.jpg

图 1.3 – 双尾假设检验中临界值与检验统计量之间的关系

根据图示，我们现在对检验统计量超过临界值表明拒绝零假设有了直观的认识。

然而，仅使用将测试统计量与假设中的临界值进行比较的方法，存在一个担忧，即测试统计量可能过大，不切实际。这种情况可能发生在结果范围很广，而这些结果被认为不属于治疗效果范围时。是否可能得到与测试统计量一样极端或更极端的结果尚不确定。为了防止错误地拒绝零假设，使用了p 值。p 值表示仅凭偶然导致观察到的值（表明拒绝零假设的值）的概率。如果 p 值相对于显著性水平较低，则可以拒绝零假设。常见的显著性水平为 0.01、0.05 和 0.10。在做出关于假设的决定之前，评估临界值与测试统计量之间的关系以及 p 值是有益的。更多内容将在第三章**，假设检验中讨论，当我们开始讨论假设检验时。

抽样策略 – 随机、系统、分层和聚类

在本节中，我们将讨论在研究中使用的不同抽样方法。总的来说，在现实世界中，由于许多原因，很难或不可能获取整个总体数据。例如，收集数据的成本在金钱和时间上都很昂贵。在许多情况下，收集所有数据是不切实际的，而且还要考虑伦理问题。从总体中抽取样本可以帮助我们克服这些问题，并且是收集数据的一种更有效的方法。通过为研究收集适当的样本，我们可以对总体属性进行统计结论或统计推断。推断性统计分析是统计思维的基本方面。本节将讨论从概率策略到在研究和工业中使用的非概率策略的不同抽样方法。

实际上存在两种类型的抽样方法：

• 概率抽样

• 非概率抽样

概率抽样

在概率抽样中，样本是根据概率理论从总体中选择的，或者是通过随机选择随机选择的。在随机选择中，总体中每个成员被选中的机会是相等的。例如，考虑一个有 10 张相似纸张的游戏。我们写上数字 1 到 10，每张纸对应一个数字。然后这些数字在一个盒子里被打乱。游戏要求随机挑选这三张纸张中的三张。因为纸张是使用相同的过程准备的，所以任何一张纸张被选中的机会（或数字 1 到 10）对每一张都是相等的。总的来说，这 10 张纸张被视为一个总体，而选出的 3 张纸张构成一个随机样本。这个例子是我们将在本章讨论的概率抽样方法之一。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_01_004.jpg

图 1.4 – 随机抽样示例

我们可以使用之前描述的（并在图 1.4**.4中展示的）抽样方法，使用numpy实现。我们将使用choice方法从给定的人口中选取三个样本。注意在choice中使用replace==False。这意味着一旦一个样本被选中，它将不会再次被考虑。注意，在以下代码中使用随机生成器是为了可重复性：

import numpy as np
# setup generator for reproducibility
random_generator = np.random.default_rng(2020)
population = np.arange(1, 10 + 1)
sample = random_generator.choice(
    population,    #sample from population
    size=3,        #number of samples to take
    replace=False  #only allow to sample individuals once
)
print(sample)
# array([1, 8, 5])

随机选择的目的在于避免当人口中的某些单位在样本中被选中的概率低于或高于其他单位时，产生偏差的结果。如今，可以通过使用计算机随机化程序来完成随机选择过程。

这里将要讨论的四种主要类型的概率抽样方法如下：

• 简单随机抽样

• 系统抽样

• 分层抽样

• 群体抽样

让我们逐一来看它们。

简单随机抽样

首先，简单随机抽样是从人口中随机选择样本的方法。通过无偏选择方法，子集（或样本）中的每个成员都有相同的机会被选中。当人口中的所有成员都与重要变量（重要特征）有相似的性质时，这种方法被使用，并且这是概率抽样的最直接方法。这种方法的优势在于最小化偏差并最大化代表性。然而，虽然这种方法有助于限制偏差的方法，但简单随机抽样存在错误的风险。这种方法也有一些局限性。例如，当人口非常大时，可能需要很高的成本和大量时间。当样本不能代表人口且研究需要再次进行此抽样过程时，需要考虑抽样误差。此外，并不是人口中的每个成员都愿意自愿参与研究，这使得获取代表大量人口的良好信息成为一个巨大的挑战。选择 10 张纸中的 3 张纸的先前例子是一个简单随机样本。

系统抽样

在这里，人口中的成员以固定的采样间隔在随机起点被选中。我们首先通过将人口中的成员数除以研究中进行的样本成员数来选择一个固定的采样间隔。然后，在采样间隔的成员数之间选择一个介于第一个成员和成员数之间的随机起点。最后，通过重复此采样过程来选择后续成员，直到收集到足够的样本。当成本和时间是研究需要考虑的主要因素时，这种方法比简单的随机抽样更快、更可取。另一方面，在简单随机抽样中，每个成员被选中的机会是均等的，而在系统抽样中，使用采样间隔规则从样本中选择一个成员进行研究。可以说，系统抽样比简单随机抽样更不随机。同样，在简单随机抽样中，成员属性与重要变量/特征的相关性也类似。让我们通过以下例子来讨论我们如何进行系统抽样。在德克萨斯州达拉斯的一所高中的班级中，有 50 名学生，但只有 10 本书要分给这些学生。通过将班级中的学生数除以书籍数（50/10 = 5）来固定采样间隔。我们还需要生成一个介于 1 和 50 之间的随机数作为随机起点。例如，取数字 18。因此，被选中获得书籍的 10 名学生如下：

18, 23, 28, 33, 38, 43, 48, 3, 8, 13

自然会提出一个问题，即区间抽样是否是一个分数。例如，如果我们有 13 本书，那么采样间隔将是 50/13 ~ 3.846。然而，我们不能选择这个分数作为代表学生数量的采样间隔。在这种情况下，我们可以选择数字 3 或 4 作为采样间隔（我们也可以选择 3 或 4 作为采样间隔）。让我们假设生成的随机起始点是 17。那么，选出的 13 名学生是这些：

17, 20, 24, 27, 31, 34, 38, 41, 45, 48, 2, 5, 9

观察前面的数字序列，在达到数字 48 之后，由于加上 4 会产生一个大于学生人数（50 名学生）的数字，因此序列从 2 重新开始（48 + 4 = 52，但由于 50 是最大值，所以我们从 2 重新开始）。因此，序列中的最后三个数字是 2、5 和 9，分别对应 4、3 和 4 的采样间隔（通过数字 50 并返回到数字 1，直到我们选出 13 名学生进行系统抽样）。

在系统抽样中，当人口成员的列表组织以匹配采样间隔时，存在偏差风险。例如，回到 50 名学生的案例，研究人员想知道学生对数学课的看法。然而，如果数学成绩最好的学生对应于数字 2、12、22、32 和 42，那么如果随机起始点是 2 且采样间隔是 10，调查可能会存在偏差。

分层抽样

这是一个基于将人口划分为同质子种群称为层的概率抽样方法。每个层根据明显不同的属性进行划分，例如性别、年龄、颜色等。这些子种群必须是不同的，以便每个层中的每个成员都有相同的机会被简单随机抽样选中。图 1.5说明了如何从两个子种群（一组数字和一组字母）中进行分层抽样以选择样本：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_01_005.jpg

图 1.5 – 分层抽样示例

以下代码示例展示了如何使用图 1.5中的示例通过numpy实现分层抽样。首先，实例被分割到相应的层：数字和字母。然后，我们使用numpy从每个层中抽取随机样本。就像之前的代码示例一样，我们利用choice方法抽取随机样本，但每个层的样本大小是基于每个层中实例的总数，而不是整个种群中的总数；例如，抽取数字的 50%和字母的 50%：

import numpy as np
# setup generator for reproducibility
random_generator = np.random.default_rng(2020)
population = [
  1, "A", 3, 4,
  5, 2, "D", 8,
  "C", 7, 6, "B"
]
# group strata
strata = {
    'number' : [],
    'string' : [],
}
for item in population:
    if isinstance(item, int):
        strata['number'].append(item)
    else:
        strata['string'].append(item)
# fraction of population to sample
sample_fraction = 0.5
# random sample from stata
sampled_strata = {}
for group in strata:
    sample_size = int(
        sample_fraction * len(strata[group])
    )
    sampled_strata[group] = random_generator.choice(
            strata[group],
            size=sample_size,
            replace=False
    )
print(sampled_strata)
#{'number': array([2, 8, 5, 1]), 'string': array(['D', 'C'], dtype='<U1')}

这种方法的主要优势是，样本中的关键群体特征更好地代表了所研究的群体，并且与整体群体成比例。这种方法有助于减少样本选择偏差。另一方面，当将群体中的每个成员分类到不同的子群体中不明显时，这种方法变得不可用。

聚类抽样

在这里，一个群体被划分为不同的子群体，称为聚类。每个聚类具有同质特征。不是在每个聚类中随机选择个别成员，而是随机选择整个聚类，并且这些聚类有同等的机会被选中作为样本的一部分。如果聚类很大，那么我们可以通过使用之前提到的抽样方法之一来选择每个聚类内的个别成员，从而进行多阶段抽样。现在讨论一个聚类抽样的例子。一家当地比萨饼店计划在其社区内扩展业务。店主想知道有多少人从他的比萨饼店订购比萨饼，以及最受欢迎的比萨饼是什么。然后他将社区划分为不同的区域，并随机选择客户形成聚类样本。他向选定的客户发送调查问卷以进行其业务研究。另一个例子与多阶段聚类抽样相关。一家连锁零售店进行一项研究，以查看连锁店中每家店铺的表现。店铺根据位置划分为子群体，然后随机选择样本形成聚类，并将样本聚类用作其店铺的表现研究。这种方法既简单又方便。然而，样本聚类并不保证能代表整个群体。

非概率抽样

另一种抽样方法是非概率抽样，在这种方法中，一个群体中的某些或所有成员没有平等的机会被选中作为样本参与研究。当无法进行随机概率抽样，且与概率抽样方法相比，这种方法更快、更容易获取数据时，会使用这种方法。使用这种方法的一个原因是因为成本和时间考虑。它允许我们通过基于便利性或某些标准进行非随机选择来轻松收集数据。这种方法可能导致比概率抽样方法更高的偏差风险。这种方法通常用于探索性和定性研究。例如，如果一组研究人员想要了解客户对其产品相关公司意见，他们会向购买并使用该产品的客户发送调查问卷。这是一种获取意见的便捷方式，但这些意见仅来自已经使用过产品的客户。因此，样本数据仅代表一组客户，不能作为公司所有客户的意见进行推广。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_01_006.png

图 1.6 – 一项调查研究的示例

之前的例子是我们想要在此处讨论的两种非概率抽样方法之一。这种方法是便利抽样。在便利抽样中，研究人员从总体中选择最易于接触的成员来形成一个样本。这种方法简单且成本低，但将获得的结果推广到整个总体是有疑问的。

配额抽样是另一种非概率抽样，其中通过非随机方式选择一个样本组，使其代表更大的总体。例如，时间有限的招聘人员可以使用配额抽样方法从专业社交网络（LinkedIn、Indeed.com 等）中寻找潜在候选人并进行面试。这种方法既经济又节省时间，但在选择过程中存在偏差。

在本节中，我们概述了概率抽样和非概率抽样。每种策略都有其优势和劣势，但它们有助于我们最小化风险，例如偏差。一个精心规划的抽样策略还将有助于减少预测模型中的错误。

摘要

在本章中，我们讨论了安装和设置 Python 环境以运行Statsmodels API 和其他必需的开源软件包。我们还讨论了总体与样本以及从样本中获得推断的要求。最后，我们解释了在统计和机器学习模型中常用的几种不同抽样方法。

在下一章中，我们将开始讨论统计分布及其对构建统计模型的影响。在第三章**假设检验中，我们将深入讨论假设检验，扩展本章观察性研究部分讨论的概念。我们还将讨论功效分析，这是一种基于现有样本数据参数和所需统计显著性水平的工具，用于确定样本量。

第二章：数据分布

在本章中，我们将涵盖数据和分布的基本方面。我们将从介绍数据类型和数据分布开始。在介绍分布的基本度量后，我们将描述正态分布及其重要特性，包括中心极限定理。最后，我们将涵盖重采样方法，如排列，以及变换方法，如对数变换。本章涵盖了开始统计建模所需的基础知识。

在本章中，我们将涵盖以下主要主题：

• 理解数据类型

• 测量和描述分布

• 正态分布和中心极限定理

• 自举

• 排列

• 变换

技术要求

本章将使用 Python 3.8。

本章的代码可以在以下位置找到 – github.com/PacktPublishing/Building-Statistical-Models-in-Python – 在ch2文件夹中。

请设置一个虚拟环境或 Anaconda 环境，并安装以下包：

• numpy==1.23.0

• scipy==1.8.1

• matplotlib==3.5.2

• pandas==1.4.2

• statsmodels==0.13.2

理解数据类型

在讨论数据分布之前，了解数据类型会有所帮助。理解数据类型至关重要，因为数据类型决定了可以使用哪些分析，因为数据类型决定了可以使用哪些操作来处理数据（这将在本章的示例中变得更加清晰）。有四种不同的数据类型：

• 名义数据

• 序数数据

• 区间数据

• 比率数据

这些类型的数据也可以分为两组。前两种数据类型（名义和序数）是定性数据，通常是非数字类别。最后两种数据类型（区间和比率）是定量数据，通常是数值。

让我们从名义数据开始。

名义数据

名义数据是对具有不同分组进行标记的数据。例如，考虑一个标志工厂中的机器。工厂通常从不同的供应商那里采购机器，这也会有不同的型号号。例如，示例工厂可能有 3 台型号 A和 5 台型号 B（见图 2.1）。机器将构成一组名义数据，其中型号 A和型号 B是不同的分组标签。对于名义数据，只有一个操作可以执行：相等。组内的每个成员都是相等的，而来自不同组的成员是不相等的。在我们的工厂示例中，一台型号 A机器将等于另一台型号 A机器，而一台型号 B机器将不等于一台型号 A机器。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_001.jpg

图 2.1 – 工厂中的两组机器

如我们所见，使用这种类型的数据，我们只能在标签下将项目分组。在下一类型的数据中，我们将介绍一个新特性：顺序。

序数数据

下一种类型的数据类似于名义数据，但表现出顺序。数据可以被标记成不同的组，并且这些组可以排序。我们称这种类型的数据为序数数据。继续以工厂为例，假设有一个Model C型号的机器，Model C与Model B由相同的供应商提供。然而，Model C是高性能版本，产生更高的输出。在这种情况下，Model B和Model C是序数数据，因为Model B是低输出机器，而Model C是高输出机器，这创造了一个自然顺序。例如，我们可以按性能的升序排列模型标签：Model B，Model C。大学教育水平也是序数数据的另一个例子，有 BS、MS 和 PhD 等级。如前所述，这种类型数据的新操作是排序，意味着数据可以被排序。因此，序数数据支持排序和相等性。虽然这种类型的数据可以按升序或降序排序，但我们不能对数据进行加减运算，这意味着Model B + Model C不是一个有意义的陈述。我们将讨论的下一类型数据将支持加减运算。

区间数据

下一种类型的数据，区间数据，用于描述存在于区间尺度上的数据，但没有明确的零点定义。这意味着两个数据点之间的差异是有意义的。以摄氏温度尺度为例。数据点是数值的，并且数据点在区间内均匀分布（例如，20 和 40 都距离 30 有 10 度）。在这个温度尺度的例子中，0 的定义是任意的。对于摄氏度来说，0 恰好被设定在水的冰点，但这是尺度设计者做出的一个任意选择。因此，区间数据类型支持相等性、排序和加减运算。

比例数据

最后一种数据类型是比例数据。与区间数据一样，比例数据是有序的数值数据，但与区间数据不同，比例数据有一个绝对零点。绝对零点意味着如果一个比例类型变量的值为零，那么该变量不存在或不存在。例如，考虑游乐园的等待时间。如果没有人在等待游乐设施，等待时间就是 0；新客人可以立即乘坐游乐设施。等待时间没有有意义的负值。0 是绝对最小值。比例数据也支持有意义的乘除运算，这使得比例数据成为支持操作最多的数据类型。

数据类型的可视化

数据可视化是理解分布和识别数据属性的关键步骤。在本章（以及整本书中），我们将使用matplotlib来可视化数据。虽然可以使用其他 Python 库进行数据可视化，但matplotlib是 Python 事实上的标准绘图库。在本节中，我们将开始使用matplotlib来可视化之前讨论的四种类型的数据。

绘制定性数据类型

由于前两种数据类型是分类的，我们将使用柱状图来可视化这些数据的分布。示例柱状图显示在图 2.2 中。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_002.jpg

图 2.2 – 柱状图中的名义数据（左侧）和有序数据（右侧）

图 2.2 左侧的柱状图显示了工厂示例中给出的模型 A机器和模型 B机器的分布。右侧的柱状图显示了工程师团队的教育水平分布示例。请注意，在教育水平柱状图中，x轴标签按教育水平从低到高排序。

生成图 2.2 所使用的代码如下所示。

代码有三个主要部分。

• 库导入：

在这个例子中，我们只从matplotlib中导入pyplot，它通常以plt的形式导入。

• 创建数据的代码：

在import语句之后，有一些语句用于创建我们将要绘制的图形所需的数据。第一个图形的数据存储在两个 Python 列表中：label和counts，分别包含机器标签和机器数量。值得注意的是，这两个列表包含相同数量的元素（两个元素）。教育数据以类似的方式存储。虽然在这个例子中，我们使用的是简单的示例数据，但在后面的章节中，我们将有额外的步骤来检索和格式化数据。

• 绘制数据的代码：

最后一步是绘制数据。由于在这个例子中我们要绘制两组数据，我们使用 subplots 方法，这将创建一个图表网格。subplots 方法的两个参数是网格中图表的行数和列数。在我们的情况下，行数是 1，列数是 2。subplots 方法返回两个对象；图，fig，和坐标轴，ax。返回的第一个对象，fig，对图有高级控制，例如保存图、在新窗口中显示图等。第二个对象，ax，将是一个单独的坐标轴对象或坐标轴对象的数组。在我们的情况下，ax 是坐标轴对象的数组——因为我们的网格有两个图表，对 ax 进行索引给我们坐标轴对象。我们使用坐标轴对象的 bar 方法创建柱状图。bar 方法有两个必需的参数。第一个必需的参数是标签列表。第二个参数是与每个标签对应的柱高，这就是为什么两个列表必须有相同长度。其他三个方法，set_title、set_ylabel 和 set_xlabel，设置相应图表属性的值：title、ylabel 和 x-label。

最后，使用 fig.show() 创建图：

import matplotlib.pyplot as plt
label = ['model A', 'model B']
counts = [3, 5]
edu_label = ['BS', 'MS', 'PhD']
edu_counts = [10, 5, 2]
fig, ax = plt.subplots(1, 2, figsize=(12, 5))
ax[0].bar(label, counts)
ax[0].set_title('Counts of Machine Models')
ax[0].set_ylabel('Count')
ax[0].set_xlabel('Machine Model')
ax[1].bar(edu_label, edu_counts)
ax[1].set_title('Counts of Education Level')
ax[1].set_ylabel('Count')
ax[1].set_xlabel('Education Level')
fig.show()

现在我们来看如何绘制其他两种数据类型的数据。

绘制定量数据类型

由于最后两种数据类型是数值型的，我们将使用直方图来可视化分布。图 2.3 中显示了两个示例直方图。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_003.jpg

图 2.3 – 柱状图中的名义数据（左侧）和有序数据（右侧）

左侧的直方图是合成等待时间数据（比率数据），可能代表游乐园的等待时间。右侧的直方图是 2022 年 4 月和 5 月达拉斯-沃斯堡地区的温度数据（来自 www.iweathernet.com/texas-dfw-weather-records)。

生成 图 2.3 所使用的代码如下。同样，代码有三个主要部分，库导入、数据创建代码和绘图代码。

与前面的例子一样，matplotlib 被导入为 plt。在这个例子中，我们还从 scipy 中导入了一个函数；然而，这个函数仅用于生成用于工作的样本数据，我们不会在这里详细讨论它。对于我们的目的，只需将 skewnorm 视为一个生成数字数组的函数即可。这个代码块与前面的代码块非常相似。

主要区别在于用于绘制数据的图表方法，hist，它创建直方图。hist方法有一个必需的参数，即要在直方图中绘制的数字序列。本例中使用的第二个参数是bins，它实际上控制了直方图的粒度 – 粒度随着更多的 bins 而增加。直方图的 bin 计数可以根据所需的视觉效果进行调整，通常通过实验设置来绘制数据：

from scipy.stats import skewnorm
import matplotlib.pyplot as plt
a = 4
x = skewnorm.rvs(a, size=3000) + 0.5
x = x[x > 0]
dfw_highs = [
    85, 87, 75, 88, 80, 86, 90, 94, 93, 92, 90, 92, 94,
    93, 97, 90, 95, 96, 96, 95, 92, 70, 79, 73, 88, 92,
    94, 93, 95, 76, 78, 86, 81, 95, 77, 71, 69, 88, 86,
    89, 84, 82, 77, 84, 81, 79, 75, 75, 91, 86, 86, 84,
    82, 68, 75, 78, 82, 83, 85]
fig, ax = plt.subplots(1,2, figsize=(12, 5))
ax[0].hist(x, bins=30)
ax[0].set_xlabel('Wait Time (hr)')
ax[0].set_ylabel('Frequency')
ax[0].set_title('Wait Times');
ax[1].hist(dfw_highs, bins=7)
ax[1].set_title('High Temperatures for DFW (4/2022-5/2022)')
ax[1].set_ylabel('Frequency')
ax[1].set_xlabel('Temperature (F)')
fig.show()

在本节中，我们简要了解了数据和分布的多样性。由于数据分布在外部世界中以许多形状和大小出现，因此拥有描述分布的方法是有用的。在下一节中，我们将讨论可用于分布的测量方法，这些测量的执行方式以及可以测量的数据类型。

测量和描述分布

在野外发现的数据分布形状和大小各异。本节将讨论如何测量分布以及哪些测量适用于四种类型的数据。这些测量将提供比较和对比不同分布的方法。本节讨论的测量可以划分为以下类别：

• 中心趋势

• 变异性

• 形状

这些测量被称为描述性统计。本节中讨论的描述性统计在数据统计摘要中常用。

测量中心趋势

中心趋势的测量有三种类型：

• 模式

• 中位数

• 均值

让我们逐一讨论它们。

模式

我们将首先讨论的中心趋势测量是模式。数据集的模式简单地说是最常见的实例。以工厂中的机器为例（见图 2.1），数据集的模式将是型号 B。在例子中，有 3 台型号 A 和 5 台型号 B，因此，使型号 B 成为最常见的 – 模式。

数据集可以是以下之一：

• 单模态 – 具有一个模式

• 多模态 – 具有多个模式

在前面的例子中，数据是单模态的。

再次使用工厂的例子，让我们想象有 3 台型号 A，5 台型号 B和 5 台型号 D（一种新型号）。那么，数据集将有两个模式：型号 B和型号 D，如图图 2.4所示。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_004.jpg

图 2.4 – 工厂中机器的多模态分布

因此，这个数据集是多模态的。

模式和数据类型

这些模式的例子使用了名义数据，但所有四种数据类型都支持模式，因为所有四种数据类型都支持等价操作。

虽然“众数”指的是最常见的实例，但在连续数据的多个峰值的多元情况下，术语“众数”通常被用来表示一个不那么严格的概念。例如，图 2.5 中的分布通常被称为多峰分布，尽管分布的峰值大小并不相同。然而，在名义数据和有序数据中，当提到分布的模态时，更常见的是使用更严格的“最常见”的定义。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/image_00_005.jpg

图 2.5 – 数据的多峰分布

现在，我们将看看如何使用 scipy 代码计算众数。scipy 库包含 stats 模块中的描述性统计函数。在这个例子中，我们从 scipy.stats 导入 mode 并计算以下数字的众数，1, 2, 3, 4, 4, 4, 5, 5：

from scipy.stats import mode
m = mode([1,2,3,4,4,4,5,5])
print(
    f"The mode is {m.mode[0]} with a count of"
    f" {m.count[0]} instances"
)
# The mode is 4 with a count of 3 instances

mode 函数返回一个包含 mode 和 count 成员的 mode 对象。不出所料，mode 和 count 成员分别包含数据集的众数和众数出现的次数。请注意，mode 和 count 成员是可索引的（就像列表一样），因为数据集可以包含多个众数。

中位数

中心测量的下一个指标是中位数。中位数是在将数值按顺序排列时出现的中间值。

中位数和数据类型

这种测量可以在有序数据、区间数据和比率数据上执行，但不能在名义数据上执行。

我们将在这里讨论两种情况。

当实例数量为奇数时寻找中位数

在 图 2.6 中显示了某些数值数据的中位数。数据被排序后，中位数被识别。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_006.jpg

图 2.6 – 使用奇数个实例识别中位数

在前面的例子中，实例数量为奇数（7 个实例），有一个中心值。然而，如果实例数量是偶数，在排序值之后，就不能简单地取中间的数字了。

当实例数量为偶数时寻找中位数

当实例数量为偶数时，取两个中间值的最小值。与众数不同，对于同一系列数据，没有多个中位数的概念。一个具有偶数个实例（8 个实例）的例子在 图 2.7 中显示。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_007.jpg

图 2.7 – 使用偶数个实例识别中位数

现在，让我们看看如何使用 numpy 计算数据集的中位数。与 scipy 一样，numpy 包含用于计算描述性统计的函数。我们将计算前面例子中列出的八个数字的中位数：

import numpy as np
values = [85, 99, 70, 71, 86, 88, 94, 105]
median = np.median(values)
print(f"The median value is {median:.2f}")
# The median value is 87.00

计算中位数的结果是 87，正如预期的那样。请注意，median 函数返回一个单一值，与前面代码示例中的 mode 函数不同。

平均值

下一个中心度量是均值，通常被称为平均值。均值由以下方程定义：

_ x  =  ∑ i=0 n  x i _ N

让我用文字解释一下这个方程。要计算均值，我们必须将所有值相加，然后将总和除以值的数量。请参考以下示例。首先将 7 个数字相加，总和达到 593。然后将这个总和除以实例数量，得到 84.7 的值。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_008.jpg

图 2.8 – 寻找均值

注意，这些值的均值和中位数（分别为 84.7 和 86）并不相同。一般来说，均值和中位数不会是相同的值，但存在一些特殊情况，均值和中位数会收敛。

均值和数据类型

至于支持的数据类型，均值适用于区间和比率数据，因为值是可以相加的。

现在，我们将看看如何使用numpy计算均值。以下代码示例显示了之前示例中值的均值的计算：

import numpy as np
values = [85, 99, 70, 71, 86, 88, 94]
mean = np.mean(values)
print(f"The mean value is {mean:.1f}")
# The mean value is 84.7

与median函数一样，mean函数返回一个单一的数字。

在结束关于中心度量的这一节之前，讨论均值和中位数在各种情况下的使用是值得的。如前所述，中位数和均值通常会是不同的值。这是由分布的形状驱动的效应。

形状对均值和中位数的影响

如果分布是对称的，均值和中位数往往会收敛。然而，如果分布不是对称的，均值和中位数会发散。

度量之间的差异程度由分布的不对称性驱动。图 2.6给出了四个示例分布来展示这一效果。分布 1 和 2 显示均值被拉向比中位数更高的值。均值会被拉向绝对值更大的值。当数据集包含（或可能包含）异常值（通常称为异常值或影响点）时，这是一个均值的重要影响，这些异常值会倾向于将其拉向它们的方向。与均值不同，如果异常值占数据比例较小，中位数不会受到异常值的影响。异常值将在测量变异性部分进一步讨论。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_009.jpg

图 2.9 – 两个不对称分布和两个对称分布

分布的下一个测量类别是变异性度量。

测量变异性

通过变异性，我们本质上是指分布的宽度。这一类别的测量如下：

• 范围

• 四分位距

• Tukey fences

• 方差

让我们讨论每一个。

范围

范围是分布中最大值和最小值之间的差值。像均值一样，范围会受到异常值的影响，因为它依赖于最大值和最小值。然而，还有一种变异方法，就像中位数一样，对异常值的存在具有鲁棒性。

让我们用代码和 numpy 来计算一个范围：

import numpy as np
values = [85, 99, 70, 71, 86, 88, 94, 105]
max_value = np.max(values)
min_value = np.min(values)
range_ = max_value - min_value
print(f"The data have a range of {range_}"
      f" with max of {max_value}"
      f" and min of {min_value}")
# The data have a range of 35 with max of 105 and min of 70

虽然 numpy 没有范围函数，但可以使用 numpy 提供的 min 和 max 函数来计算范围。

四分位距

下一个变异度量是通过排序数据然后将数据分成四个相等的部分来确定的。四个部分的边界是四分位数，它们被称为以下：

• 下四分位数（Q1）

• 中间四分位数（Q2）

• 上四分位数（Q3）

以下是一个四分位数的例子。像中位数一样，只要异常值占数据集的比例很小，四分位数对异常值就具有鲁棒性。请注意，中间四分位数实际上就是中位数。一种对异常值比 范围 部分中讨论的正常范围更不敏感的调整范围测量方法是中间四分位数，即 四分位距（IQR）。IQR 是上四分位数（Q3）和下四分位数（Q1）之间的差值。虽然这个范围对异常值不太敏感，但它只包含 50% 的数据。因此，四分位距很可能是不太能代表数据总变异的。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_010.jpg

图 2.10 – Q1, Q2 和 Q3

我们可以使用 numpy 和 scipy 来计算四分位数和 IQR 范围。在下面的代码示例中，我们使用 quantiles 函数来计算四分位数。我们在这里不会讨论 quantiles，除了提到 quantiles 是一种泛化，其中数据可以被分成任何数量的相等部分。由于我们是将数据分成四个相等的部分来计算四分位数，所以用于计算的 quantiles 值是 0.25、0.5 和 0.75。然后可以使用 Q1 和 Q3 来计算 IQR。但是，我们也可以使用 scipy 中的 iqr 函数来进行计算：

import numpy as np
from scipy import stats
values = [85, 99, 70, 71, 86, 88, 94]
quartiles = np.quantile(values, [0.25, 0.5, 0.75],
    method="closest_observation")
print(f"The quartiles are Q1: {quartiles[0]},
    Q2: {quartiles[1]}, Q3: {quartiles[2]}")
iqr = stats.iqr(values,interpolation='closest_observation')
print(f"The interquartile range is {iqr}")
# The quartiles are Q1: 71, Q2: 85, Q3: 88
# The interquartile range is 17

注意 quantiles 函数和 iqr 函数中 method 和 interpolation 关键字参数的使用。可以对这些关键字参数使用几种选项，这将导致不同的结果。

四分位数通常用箱线图来可视化。以下 图 2*.11 展示了箱线图的主要部分。箱线图由两个主要部分组成：

• 箱体

• 须线

箱体部分代表包含在 IQR 内的 50% 的数据。须线从箱体的边缘延伸到长度为 k * IQR，其中 k 通常选择为 1.5。任何超出须线的值都被认为是异常值。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_011.jpg

图 2.11 – 箱线和须图的部分

图 2.12展示了直方图和箱线图如何可视化对称和非对称分布的变异。注意非对称数据的箱线图在左侧被压缩，在右侧被扩展，而另一个箱线图则明显对称。虽然箱线图有助于可视化数据的对称性和异常值的存在，但分布的模态可能不会明显。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_012.jpg

图 2.12 – 对称和非对称分布的箱线图和直方图比较

在探索时，通常使用多种可视化方式，因为每种类型的可视化都有其自身的优势和劣势。通常使用多种可视化方式，因为每种类型的可视化都有其自身的优势和劣势。

Tukey fences

在最后几节关于测量的内容中，异常值的概念出现了几次。异常值是与主要分布相比不典型或异常的值。虽然有一些方法可以将数据点分类为异常值，但没有一种通用的方法可以用来分类数据点为异常值。定义异常值通常应该根据数据分析的使用情况进行，因为根据应用领域会有不同的因素需要考虑。然而，值得提及的是，在箱线图示例中展示的常见技术，称为 Tukey fences。下限和上限 Tukey fences 基于四分位数间距（IQR）并定义如下：

• 下限围栏：Q1 − k(IQR)

• 上限围栏：Q3 + k(IQR)

如前所述，k 通常被选为默认值 1.5，但针对特定应用领域可能有一个更合适的值。

现在我们来看一下如何使用numpy和scipy计算 Tukey fences。这个代码示例将基于前面的示例，因为没有直接计算围栏的函数。我们再次使用numpy和scipy计算四分位数和四分位数间距（IQR）。然后，我们将这些操作应用于前面方程中列出的值：

import numpy as np
from scipy import stats
values = stats.norm.rvs(10, size=3000)
q1, q3 = np.quantile(values, [.25, .75],
    method='closest_observation')
iqr = stats.iqr(values,interpolation='closest_observation')
lower_fence = q1 - iqr * 1.5
upper_fence = q3 + iqr * 1.5
# may vary due to randomness in data generation
print(f"The lower fence is {lower_fence:.2f} and the upper
    fence is {upper_fence:.2f}")
# The lower fence is 7.36 and the upper fence is 12.67

在这个例子中，我们使用了numpy和scipy；然而，如前所述，scipy的计算可以用Q3-Q1来替换。

方差

在本节中将要介绍的最后一个变异度量是方差。方差是衡量离散程度的度量，可以理解为数值与平均值之间的分散程度。方差的公式，用 S 2 表示，如下所示：

S 2 =  ∑ (x i −  _ x ) 2 _ N − 1

在这个方程中，项（x i −  _ x ）被认为是与平均值的偏差，这导致另一个与方差密切相关的度量——标准差，它是方差的平方根。标准差的公式，用σ表示，如下给出：

σ = √ _ S 2  = √ ___________  ∑ (x i −  _ x ) 2 _ N − 1

通常情况下，分布范围越宽，其方差和标准差也会越大，但这些值不如范围或四分位距（IQR）那样容易解释。这些概念将在下一节中更详细地介绍，在正态分布的背景下，这将提供对这些值所测量的更清晰的直观理解。

同样，这些值将通过代码使用numpy进行计算。方差和标准差的函数分别是var和std：

import numpy as np
values = [85, 99, 70, 71, 86, 88, 94]
variance = np.var(values)
standard_dev = np.std(values)
print(f"The variance is {variance:.2f} and the standard
    deviation is {standard_dev:.2f}")
# The variance is 101.06 and the standard deviation is 10.05

测量形状

下一种类型的测量与分布的形状有关。具体如下：

• 偏度

• 偏度

让我们逐一讨论它们。

偏度

第一项测量是偏度。简单来说，偏度是测量不对称性[1]。图 2.13展示了偏斜分布的例子。

存在两种类型的偏斜分布：

• 左偏斜

• 右偏斜

分布在主导尾的方向上偏斜，意味着具有右侧主导尾的分布是右偏斜的，而具有左侧主导尾的分布是左偏斜的（如图 2.13所示）。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_013.jpg

图 2.13 – 展示偏度的分布

由于现代软件包可以轻易计算，因此此处不会展示偏度的公式。偏度计算的结果可以用来确定偏度和偏斜的方向。如果偏度值为 0 或接近 0，则分布没有表现出强烈的偏斜。如果偏度值为正，则分布是右偏斜的；如果偏度值为负，则分布是左偏斜的。偏度值的绝对值越大，分布的偏斜程度越明显。以下代码示例展示了如何使用scipy计算偏度：

from scipy.stats import skewnorm, norm
from scipy.stats import skew as skew_calc
# generate data
skew_left = -skewnorm.rvs(10, size=3000) + 4
skew_right = skewnorm.rvs(10, size=3000) + 3
symmetric = norm.rvs(10, size=3000)
# calculate skewness
skew_left_value = skew_calc(skew_left)
skew_right_value = skew_calc(skew_right)
symmetric_value = skew_calc(symmetric)
# Output may vary some due to randomness of generated data
print(f"The skewness value of this left skewed
    distribution is {skew_left_value:.3f}")
print(f"The skewness value of this right skewed
    distribution is {skew_right_value:.3f}")
print(f"The skewness value of this symmetric distribution
    is {symmetric_value:.3f}")

本节中涵盖的另一种形状测量是峰度。

偏度

峰度是衡量分布尾部相对于正态分布的轻重程度[2]。虽然正态分布尚未深入探讨，但峰度的概念仍然可以讨论。轻尾分布意味着更多的数据接近或围绕分布的众数。相反，重尾分布意味着比接近众数的数据，更多的数据位于分布的边缘。图 2.14展示了轻尾分布、正态分布和重尾分布。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/image_00_014.jpg

图 2.14 – 以正态分布为参考的尾部分布

由于现代软件包可以轻易地计算出来，这里将不展示峰度的公式。如果峰度值为 0 或接近 0，则分布不显示峰度。如果峰度值为负，则分布显示轻尾性；如果峰度值为正，则分布显示重尾性。以下代码示例展示了如何使用scipy计算峰度：

from scipy.stats import norm
from scipy.stats import gennorm
from scipy.stats import kurtosis
# generate data
light_tailed = gennorm.rvs(5, size=3000)
symmetric = norm.rvs(10, size=3000)
heavy_tailed = gennorm.rvs(1, size=3000)
# calculate skewness
light_tailed_value = kurtosis(light_tailed)
heavy_tailed_value = kurtosis(heavy_tailed)
symmetric_value = kurtosis(symmetric)
# Output may vary some due to randomness of generated data
print(f"The kurtosis value of this light-tailed
    distribution is {light_tailed_value:.3f}")
print(f"The kurtosis value of this heavy_tailed
    distribution is {heavy_tailed_value:.3f}")
print(f"The kurtosis value of this normal
    distribution is {symmetric_value:.3f}")

在本节中，我们介绍了用于测量和描述数据分布的常见描述性统计方法。这些测量提供了描述和比较分布的共同语言。本章讨论的概念是未来章节中许多概念的基础。在下一节中，我们将讨论正态分布，并使用这些测量来描述正态分布。

正态分布和中心极限定理

讨论正态分布时，我们指的是钟形的标准正态分布，它正式同义于高斯分布，以 18 世纪和 19 世纪的数学家和物理学家卡尔·弗里德里希·高斯的名字命名。高斯在许多方面做出了贡献，包括近似概念，并在 1795 年发明了最小二乘法和正态分布，这些在统计建模技术（如最小二乘回归）中常用[3]。标准正态分布也被称为参数分布，其特征是分布对称，数据点分散的概率在均值周围一致——也就是说，数据出现在均值附近比出现在远离均值的频率更高。由于这个分布中的位置数据遵循概率定律，我们可以称之为标准正态概率分布。顺便提一下，在统计学中，不是概率分布的分布是通过基于非随机选择的非概率抽样生成的，而概率分布是基于随机抽样生成的。基于概率和非概率的分布都可以有标准正态分布。标准正态分布既不偏斜也不显示峰度。它在整个分布中具有相同的方差，并且在自然界中经常出现。经验法则用于描述这个分布具有围绕均值μ的三个相关的标准差。关于这个分布有两个不同的假设：

• 第一、第二和第三标准差分别包含 68%、95%和 99.7%的测量值。

• 均值、中位数和众数都相等

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_015.jpg

图 2.15 – 标准正态分布

正态分布的两种常见形式如下：

• 概率密度分布

• 累积密度分布

如前所述，概率密度分布基于随机抽样，而累积密度分布基于累积数据，这些数据不一定随机。

标准正态分布的双尾概率密度函数如下：

f(x) =  e −(x−μ) 2 _ 2 σ 2  _ σ √ _ 2π

标准正态分布的左尾累积函数如下：

f(x) =  ∫ −∞ x e −x 2 _  2 _ √ _ 2π

在统计建模方面，正态分布代表着平衡和对称。这在构建统计模型时非常重要，因为许多模型都假设正态分布，并且对许多偏离该假设的情况不稳健，因为它们是围绕平均值构建的。因此，如果此类模型中的变量不是正态分布的，模型的误差将会增加且不一致，从而降低模型稳定性。当考虑统计模型中的多个变量时，如果两者都是正态分布的，它们的相互作用更容易被近似。

在以下 图 2*.16* 中，左边的图中，变量 X 和 Y 相互作用并在平均值周围形成集中分散。在这种情况下，使用 X 和均值线或线性距离对 Y 进行建模可以做得相当好。然而，如果两个变量的分布是偏斜的，如右边的图所示，这会导致两个变量之间的方差不是常数，从而导致误差分布不均，输出不可靠。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_16.jpg

图 2.16 – 双变量正态分布（左）和偏斜分布（右）

在线性分类和回归模型的情况下，这意味着一些结果会比其他结果更好，而一些结果可能会非常糟糕。有时使用基本的模型指标很难评估这一点，需要更深入地分析模型以防止相信可能最终会误导的结果。此外，部署到生产环境将非常危险。关于这一点将在 第六章 中进一步讨论。

中心极限定理

在采样数据时，经常会遇到非正态数据的问题。这可能是由于多个原因造成的，例如总体没有正态分布或者样本不能代表总体。中心极限定理，在统计推断中非常重要，它假设如果从具有特定均值μ和标准差σ的总体中随机抽取 n 个观察值的样本，那么从随机选择的子样本分布的均值构建的抽样分布将接近一个具有大致相同的均值μ和标准差σ的正态分布，计算公式为√ _ ∑ (x i − μ) 2 _ N  ，与总体相同。下一节将使用自助法来展示中心极限定理的实际应用。在后面的章节中，将讨论变换技术，提供重塑不符合正态分布的数据分布的方法，以便仍然可以有效地应用需要正态分布的工具。

自助法

自助法是一种重采样方法，它使用随机抽样——通常是放回抽样——通过从样本分布的子集中重采样来生成关于总体的统计估计，如下所示：

• 置信区间

• 标准误

• 相关系数（皮尔逊相关系数）

理念是，通过重复采样样本分布的不同随机子集并每次取平均值，给定足够的重复次数，将开始使用每个子样本的平均值来近似真实的总体。这直接遵循中心极限定理的概念，即重申，采样均值开始近似围绕原始分布均值的正态采样分布，随着样本量和计数增加。当分布中样本数量相对于特定测试所需的样本量有限时，自助法很有用，但需要进行推断。

如第一章中所述，采样与泛化，时间和费用等限制是获取样本而不是总体的常见原因。因为自助法的基本概念是使用样本对总体进行假设，所以当总体的真实统计参数（如百分位数和方差）已知时，将此技术应用于总体并不有利。至于样本制备，样本中属性的平衡应代表总体的真实近似。否则，结果可能会误导。例如，如果动物园中物种的总体是 40%的爬行动物和 60%的哺乳动物，而我们想通过自助法来估计它们的寿命并确定其置信区间，那么必须确保应用于自助法的数据集包含 40%的爬行动物和 60%的哺乳动物的分布；例如，15%的爬行动物和 85%的哺乳动物的分布会导致误导性结果。换句话说，样本分层应与总体的比例相平衡。

置信区间

如前所述，自助法的一个有用应用是围绕稀疏定义或有限的数据集创建置信区间——也就是说，具有广泛值范围但没有许多样本的数据集。考虑一个使用statsmodels中的"Duncan"数据集进行假设检验的例子，该数据集包含按职业、类型、教育和声望划分的收入。虽然这是一个完整的数据集，但考虑到采样方法未提及，并且不太可能考虑所有职业和类型的所有工人的所有收入，我们可以将此数据集视为样本。为了获取数据集，我们首先加载matplotlib、statsmodels、pandas和numpy库。然后我们下载数据集并将其存储在pandas DataFrame 的df_duncan变量中。在此之后，我们将“prof"、“wc"和“bc"类型分别重新编码为“professional"、“white-collar"和“blue collar"。最后，我们创建了两个独立的pandas DataFrame；一个用于专业工作类型，另一个用于蓝领工作类型，因为这是我们打算使用自助法分析的两组子集：

import matplotlib.pyplot as plt, statsmodels.api as sm, pandas as pd, numpy as np, scipy.stats
df_duncan = sm.datasets.get_rdataset("Duncan",
    "carData").data
df_duncan.loc[df_duncan['type'] == 'prof',
    'type'] = 'professional'
df_duncan.loc[df_duncan['type'] == 'wc',
    'type'] = 'white-collar'
df_duncan.loc[df_duncan['type'] == 'bc',
    'type'] = 'blue-collar'
df_professional = df_duncan.loc[(
    df_duncan['type'] == 'professional')]
df_blue_collar = df_duncan.loc[(
    df_duncan['type'] == 'blue-collar')]

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_017.jpg

图 2.17 – 显示 statsmodels Duncan 数据的前五行表格

然后，我们构建了一组绘图函数，如接下来所示。在 plot_distributions() 中，我们表示 p=5，这意味着在显著性水平为 0.05（1.00 - 0.05 = 0.95，因此，95% 置信度）时，p 值将是显著的。然后我们将此值除以 2，因为这将是一个双尾测试，意味着我们想要知道整个区间而不是仅仅一个边界（在 第一章 中作为代表性统计量的讨论）。至于绘图，我们使用 matplotlib 中的直方图（hist() 函数）可视化数据，然后使用我们用 numpy 函数 percentile() 构建的 axvline() 函数绘制 95% 的抽样置信区间。

自助抽样中的百分位数

当应用于原始数据时，百分位数只是那样，但当应用于自助抽样分布时，它就是置信区间。

要简单地说明置信区间，95% 置信区间意味着在取出的每个 100 个样本均值中，有 95 个将落在该区间内。在 numpy 的 percentile() 函数中，我们使用 p=5 来支持 1-p 是置信水平，其中 p 是显著性水平（想想 p-value，其中任何等于或低于 p 的值都是显著的）。由于测试是双尾的，我们将 p 除以 2，并在左尾和右尾各分割 2.5，因为我们有一个对称的标准正态分布。subplot(2,1,...) 代码创建了两行一列。图形的轴 0 用于专业收入，轴 1 用于蓝领收入：

def plot_distributions(n_replicas, professional_sample, blue_collar_sample, professional_label, blue_collar_label, p=5):
    fig, ax = plt.subplots(2, 1, figsize=(10,8))
    ax[0].hist(professional_sample, alpha=.3, bins=20)
    ax[0].axvline(professional_sample.mean(),
        color='black', linewidth=5)
# sampling distribution mean
    ax[0].axvline(np.percentile(professional_sample, p/2.),
        color='red', linewidth=3, alpha=0.99)
# 95% CI Lower limit (if bootstrapping)
    ax[0].axvline(np.percentile(professional_sample,
        100-p/2.), color='red', linewidth=3, alpha=0.99)
# 95% CI Upper Limit  (if bootstrapping)
    ax[0].title.set_text(str(professional_label) +
        "\nn = {} Resamples".format(n_replicas))
    ax[1].hist(blue_collar_sample, alpha=.3, bins=20)
    ax[1].axvline(blue_collar_sample.mean(), color='black',
        linewidth=5) # sampling distribution mean
    ax[1].axvline(np.percentile(blue_collar_sample, p/2.),
        color='red', linewidth=3, alpha=0.99)
# 95% CI Lower limit (if bootstrapping)
    ax[1].axvline(np.percentile(blue_collar_sample,
        100-p/2.), color='red', linewidth=3, alpha=0.99)
# 95% CI Upper Limit (if bootstrapping)
    ax[1].title.set_text(str(blue_collar_label) +
        "\nn = {} Resamples".format(n_replicas))
    if n_replicas > 1:
        print("Lower confidence interval limit: ",
            np.percentile(round(professional_sample,4),
            p/2.))
        print("Upper confidence interval limit: ",
            np.percentile(round(professional_sample,4),
            100-p/2.))
        print("Mean: ", round(professional_sample,
            4).mean())
        print("Standard Error: ",
            round(professional_sample.std() /
            np.sqrt(n_replicas), 4) )
        print("Lower confidence interval limit: ",
            np.percentile(round(blue_collar_sample,4),
            p/2.))
        print("Upper confidence interval limit: ",
            np.percentile(round(blue_collar_sample,4),
            100-p/2.))
        print("Mean: ", round(blue_collar_sample,4).mean())
        print("Standard Error: ",
            round(blue_collar_sample.std() /
            np.sqrt(n_replicas), 4) )
    else:
        print("At least two samples required to create the following statistics:\nConfidence Intervals\nMean\nStandard Error")

在原始数据集中，有 18 个 professional 职业类型的收入数据点和 21 个 blue-collar 职业类型的收入数据点。专业职业类型的 95% 置信区间从 29.50 到 79.15，平均值为 60.06。该区间对于蓝领职业类型从 7.00 到 64.00，平均值为 23.76。根据 图 2.18，收入差异之间存在合理的重叠，这导致了重叠的置信区间。因此，可以合理地假设蓝领和专业人士之间的收入没有统计学上的显著差异。然而，这个数据集的样本量非常有限：

n_replicas=0
plot_distributions(n_replicas=n_replicas,
professional_sample=df_professional['income'],
    blue_collar_sample=df_blue_collar['income'],
    professional_label="Professional",
    blue_collar_label="Blue Collar")

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_018.jpg

图 2.18 – 带有 95% 百分位数线的原始数据分布

在以下代码中，使用 pandas 的 .sample() 函数，我们随机重新抽取每个分布中 50%（frac=0.5）的收入值 1,000 次，并每次计算一个新的均值，将其附加到以 _bootstrap_means 结尾的 Python 列表中。使用这些列表，我们推导出新的 95%置信区间。图 2*.19* 显示，与数据集中的标准差和收入值相比，使用每个重新抽样的子集的平均值的新样本分布。replace=True 参数允许多次重新抽样相同的记录（在随机发生的情况下），这是 bootstrap 的要求。

在执行 bootstrap 过程之后，我们可以看到收入已经开始以大致标准正态、高斯形式分布。值得注意的是，从这个实验中，置信区间不再重叠。专业和蓝领群体之间置信区间分离的含义是，在 95%的置信水平下，可以证明两种工作类型的收入之间存在统计学上的显著差异。专业收入水平的置信区间现在是 48.66 到 69.89，平均值为 60.04，而蓝领的置信区间是 14.60 到 35.90，平均值为 23.69：

n_replicas = 1000
professional_bootstrap_means = pd.Series(
    [df_professional.sample(frac=0.5, replace=True)
    ['income'].mean() for i in range(n_replicas)])
blue_collar_bootstrap_means = pd.Series(
    [df_blue_collar.sample(frac=0.5, replace=True)
    ['income'].mean() for i in range(n_replicas)])

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_019.jpg

图 2.19 – 1,000 次 bootstrap 抽样的 1,000 个样本均值的 95%置信区间分布

在这里，你可以注意到分布更紧密地聚集在平均值周围，置信区间更紧密。

如前所述，bootstrap 可以用来获得分布的除置信区间之外的不同统计参数。

标准误差

另一个常用的指标是标准误差，σ √ n。我们可以使用最后几个变量，professional_bootstrap_means 和 blue\_collar\_bootstrap\_means 来计算这个值，因为这些变量包含了通过 bootstrap 过程获得的新均值分布。我们还可以看到，标准误差——通过将标准差除以样本数量的平方根（在我们的情况下，n\_replicas，代表从每个随机重新抽样的样本中获得的平均值数量）——随着重新抽样的体积增加而减小。我们使用以下代码来计算专业和蓝领类型收入 bootstrap 均值的标准误差。下表，图 2*.20* 显示，随着 n 的增加，标准误差降低：

scipy.stats.sem(professional_bootstrap_means)
scipy.stats.sem(blue_collar_bootstrap_means)

n	专业 标准误差	蓝领 标准误差
10 个副本	0.93	2.09
10,000 个副本	0.03	0.04

图 2.20 – n = 10 和 n = 10,000 bootstrap 副本的标准误差表

Bootstrap 的另一个用例是皮尔逊相关系数，我们将在下一节讨论。

相关系数（皮尔逊相关）

通常，使用小样本量很难找到，因为相关性取决于两个变量的协方差。随着变量重叠得更加显著，它们的相关性就越高。然而，如果重叠是样本量小或抽样误差的结果，这种相关性可能是有代表性的。图 2.21显示了不同引导子样本计数的相关性表。随着分布形成更多的原生区分，相关性从小的正相关减少到接近零的量。

要测试来自原始数据集的 10 个记录样本的相关性，请参阅以下内容：

df_prof_corr = df_professional.sample(n=10)
df_blue_corr = df_blue_collar.sample(n=10)
corr, _ = scipy.stats.pearsonr(df_prof_corr['income'],
    df_blue_corr['income'])

要测试引导平均值样本的相关性：

n_replicas = n_replicas
professional_bootstrap_means = pd.Series([df_prof_corr.sample(frac=0.5,replace=False).income.mean()for i in range(n_replicas)])
blue_collar_bootstrap_means = pd.Series([df_blue_corr.sample(frac=0.5, replace=False).income.mean() for i in range(n_replicas)])
corr, _ = scipy.stats.pearsonr(
    professional_bootstrap_means,
    blue_collar_bootstrap_means)
print(corr)

n	Pearson 相关系数
原始数据中的 10 个样本	0.32
10 个副本	0.22
10,000 个副本	-0.003

图 2.21 – Pearson 相关系数表与原始样本并排

通常，运行大约 1,000 到 10,000 个引导副本是很常见的。然而，这取决于正在引导的数据类型。例如，如果从人类基因组序列数据集中引导数据，引导样本 1 亿次可能是有用的，但如果引导一个简单的数据集，引导 1,000 次或更少可能是有用的。最终，研究人员应进行视觉检查平均值的分布，以确定结果与预期相比是否合理。与统计学一样，最好有一些领域知识或专业知识来帮助验证发现，因为这可能是决定引导复制次数的最佳方法。

引导程序（Bootstrapping）也应用于机器学习，它是自助聚合（bootstrap aggregation）概念的基础，也称为Bagging，这是一个结合基于自助子样本分布构建的预测模型输出的过程。随机森林（Random Forest）是执行此操作的一种流行算法。在 Bagging 算法中，引导程序的目的在于保持非参数（将在后续章节中进一步讨论）分类的低偏差行为，同时减少方差，因此将引导程序作为一种最小化建模误差中偏差-方差权衡重要性的方法。

在下一节中，我们将考虑另一种非参数测试，即使用重采样数据的排列测试。

排列

在跳入这个测试分析之前，我们将回顾一些排列组合的基本知识。

排列和组合

排列和组合是两种数学技术，用于从一组对象中创建子集，但以两种不同的方式。在排列中，对象的顺序很重要，而在组合中则不重要。

为了轻松理解这些概念，我们将考虑两个例子。晚会上有 10 个人。派对的组织者想要随机给 3 个人颁发价值 1000 美元、500 美元和 200 美元的奖品。问题是*有多少种分配奖品的方式？*另一个例子是，组织者将从派对的 10 个人中选出 3 个人，每人颁发价值 500 美元的等额奖品。组织者并不关心这 3 个被选中的人中谁得到哪个奖品。在两个例子中，胡伊、保罗和斯图尔特是我们的获奖者，但在第一个例子中，可能会有不同的结果，例如，如果保罗赢得 200 美元的奖品、500 美元的奖品或 1000 美元的奖品。

**$**1,000	**$**500	**$**200
胡伊	保罗	斯图尔特
保罗	胡伊	斯图尔特
保罗	斯图尔特	胡伊
胡伊	斯图尔特	保罗
斯图尔特	胡伊	保罗
斯图尔特	保罗	胡伊

图 2.22 – 分配给胡伊、保罗和斯图尔特的奖品分布表

然而，在第二个例子中，由于 3 个奖品的价值相同，为 500 美元，奖品排列的顺序并不重要。

让我们更仔细地看看这两个排列组合的例子。第一个例子是一个排列例子。由于有 10 个人，我们有 10 种可能性从这 10 个人中选择一个人来颁发 1000 美元的奖品。如果这个人被选中赢得 1000 美元的奖品，那么选择另一个人颁发 500 美元的奖品就只剩下 9 种可能性，最后，我们有 8 种可能性从这 10 个人中选择一个人来颁发 200 美元的奖品。然后，我们有 1098 = 720 种分配奖品的方式。排列的数学公式如下：

P(n, r) = n! / (n-r)!

在这里，P(n, r)是排列的数量，n 是集合中的总对象数，r 是从集合中选择的对象数。在这个例子中，n = 10 且 r = 3，所以我们看到：

P(10,3) = 10! / (10-3)! = 109 * 87 * 65 * 43 * 21 / 76 * 54 * 32 * 1 = 10*9 * 8 = 720

从晚会上 10 个人中选出 3 个人来分配价值 1000 美元、500 美元和 200 美元的 3 个奖品，共有 720 种选择方式。

在 Python 中，有一个名为itertools的包可以帮助我们直接找到排列。读者可以查看以下链接 – docs.python.org/3/library/itertools.xhtml – 获取有关此包的更多信息。我们需要将此包导入 Python 环境以进行排列：

from itertools import permutations
# list of 10 people in the party
people = ['P1','P2','P3','P4','P5','P6','P7','P8','P9','P10']
# all the ways that the 3 prizes are distributed
perm = permutations(people, 3)
list_perm = list(perm)
print(f"There are {len(list_perm)} ways to distribute the prizes!")

在前面的 Python 代码中，我们创建了一个包含 10 个人的列表 people，分别是 P1 到 P10，然后使用 itertools 中的 permutations 函数来获取所有分配奖项的方式。这种方法接受一个包含 10 个人的列表作为输入，并返回一个包含所有可能性的元组列表，即从这 10 个人中选择 3 个人来分配 1,000 美元、500 美元和 200 美元的奖项。因为共有 720 种分配奖项的方式，这里我们只打印 Python 代码生成的 10 种第一种方式：

print(f"The 10 first ways to distribute the prizes: \n
    {list_perm[:10]} ")

前面代码的输出是分配奖项的 10 种第一种方式：

[('P1', 'P2', 'P3'), ('P1', 'P2', 'P4'), ('P1', 'P2', 'P5'), ('P1', 'P2', 'P6'), ('P1', '``P2', 'P7')]

如果我们有 10 件不同的礼物，每个参加派对的人都可以带一件礼物回家。有多少种分配这些礼物的方式？共有 3,628,800 种方式。这是一个非常大的数字！读者可以用以下代码进行验证：

#list of 10 people in the party
people = ['P1','P2','P3','P4','P5','P6','P7','P8','P9','P10']
# all the ways that the 10 different gifts are distributed
perm = permutations(people)
list_perm = list(perm)
print(f"There are {len(list_perm)}
    ways to distributed the gifts!")

回到第二个例子，因为 3 个奖项的价值都是 500 美元，所以 3 个被选中的人的顺序并不重要。那么，如果选中的 3 个人是 Huy、Paul 和 Stuart，如 图 2*.22* 所示，在第一个例子中分配奖项的方式有 6 种。然后，将相同数量的 500 美元分配给 Huy、Paul 和 Stuart 的方式只有 1 种。组合的数学公式如下：

C(n, r) = n! / (r! * (n - r)!)

在这里，C(n, r) 表示组合数，n 是集合中对象的总数，r 是可以从集合中选择的对象数量。同样，我们可以计算出有

10! / (3! * (10 - 3)!) = 10 * 9 * 8 / (1 * 2 * 3) = 720 / 6 = 120

分配 3 个价值 500 美元的奖项的方式。

在 Python 中，我们同样使用 itertools 包，但不是使用 permutations 函数，而是导入 combinations 函数：

 from itertools import combinations
# list of 10 people in the party
people = ['P1','P2','P3','P4','P5','P6','P7','P8','P9','P10']
# all the ways that the 3 prizes are distributed
comb = combinations(people, 3)
list_comb = list(comb)
print(f"There are {len(list_comb)} ways to distribute the prizes!")

排列测试

排列测试是一种非参数测试，它不要求数据呈正态分布。自助法和排列法都适用于重采样技术，但适用于不同的用途，一个用于估计统计参数（自助法），另一个用于假设检验。排列测试用于测试来自同一总体生成的两个样本之间的零假设。它有不同的名称，如 精确测试、随机化测试和重新随机化测试。

在实现 Python 代码之前，我们先来看一个简单的例子，以便更好地理解。我们假设有两组人，一组代表儿童（A），另一组代表 40 岁以上的人（B），如下所示：

A = [3,5,4] 和 B = [43,41,56,78,54]

样本 A 和 B 之间年龄的平均差异是

43 + 41 + 56 + 78 + 54 __________________ 5 - 3 + 5 + 4 * 3 = 50.4

我们将 A 和 B 合并成一个单一集合，记为 P，如下所示：

P = [3,5, 4,43,41,56,78,54].

然后，我们取 P 的一个排列，例如以下内容：

P_new = [3,54, 78, 41, 4, 43, 5, 56]

然后，我们将P_new重新划分为两个子集，分别称为A_new和B_new，它们的大小分别与 A 和 B 相同：

A_new = [3,54,78] 和 B_new = [41,4,43,5,56]

然后，A_new和B_new之间的年龄均值差异为 15.2，这低于 A 和 B 之间的原始年龄均值差异（50.4）。换句话说，排列后的P_new对 p 值没有贡献。我们可以观察到，从 P 的所有可能排列中抽取的只有一个排列的均值差异大于或等于原始的均值差异 P。现在我们将代码实现于 Python 中：

import numpy as np
# create permutation testing function
def permutation_testing(A,B,n_iter=1000):
#A, B are 2 lists of samples to test the hypothesis,
#n_iter is number of iterations with the default is 1000
    differences = []
    P = np.array(A+B)
    original_mean = np.array(A).mean()- np.array(B).mean()
    for i in range(n_iter):
      np.random.shuffle(P)#create a random permutation of P
      A_new = P[:len(A)] # having the same size of A
      B_new = P[-len(B):] # having the same size of B
      differences.append(A_new.mean()-B_new.mean())
    #Calculate p_value
    p_value = round(1-(float(len(np.where(
        differences<=original_mean)[0]))/float(n_iter)),2)
    return p_value

在前面的 Python 代码中，A 和 B 是两个样本，我们想知道它们是否来自同一个更大的总体；n_ter是我们想要执行的迭代次数；这里，1,000 是默认的迭代次数。

让我们以 10,000 次迭代对示例中的两组人员进行置换检验：

A = [3,5,4]
B = [43,41,56,78,54]
permutation_testing(A,B,n_iter=10000)

得到的 p 值为 0.98。这意味着我们未能拒绝零假设，或者说没有足够的证据来确认样本 A 和 B 来自同一个更大的总体。

接下来，我们将探讨许多需要正态分布假设的统计测试中的一个重要且必要的步骤。

转换

在本节中，我们将考虑三种转换：

• 对数转换

• 平方根转换

• 立方根转换

首先，我们将导入numpy包以创建一个从 Beta 分布中抽取的随机样本。Beta 分布的文档可以在这里找到：

numpy.org/doc/stable/reference/random/generated/numpy.random.beta.xhtml

样本df有 10,000 个值。我们还使用matplotlib.pyplot创建不同的直方图。其次，我们通过使用对数转换、平方根转换和立方根转换来转换原始数据，并绘制四个直方图：

import numpy as np
import matplotlib.pyplot as plt
np.random.seed(42) # for reproducible purpose
# create a random data
df = np.random.beta(a=1, b=10, size = 10000)
df_log = np.log(df) #log transformation
df_sqrt = np.sqrt(df) # Square Root transformation
df_cbrt = np.cbrt(df) # Cube Root transformation
plt.figure(figsize = (10,10))
plt.subplot(2,2,1)
plt.hist(df)
plt.title("Original Data")
plt.subplot(2,2,2)
plt.hist(df_log)
plt.title("Log Transformation")
plt.subplot(2,2,3)
plt.hist(df_sqrt)
plt.title("Square Root Transformation")
plt.subplot(2,2,4)
plt.hist(df_cbrt)
plt.title("Cube Root Transformation")
plt.show()

以下为代码的输出：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_02_023.jpg

图 2.23 – 原始数据和转换数据的直方图

通过转换，我们可以看到转换后的直方图比原始直方图更接近正态分布。似乎在这个例子中，最佳转换是立方根转换。在现实世界的数据中，确定是否需要转换以及应该使用哪种转换非常重要。

其他数据转换方法，例如查找重复数据、处理缺失值和特征缩放，将在以下章节的 Python 实际应用案例中讨论。

摘要

在本章的第一节中，我们学习了数据的类型以及如何可视化这些类型的数据。然后，我们介绍了如何描述和测量数据分布的属性。我们学习了标准正态分布，为什么它很重要，以及通过展示自助法来展示中心极限定理在实践中的应用。我们还学习了如何利用自助法利用非正态分布的数据，通过置信区间来测试假设。接下来，我们介绍了排列组合等数学知识，并介绍了排列测试作为自助法之外的另一种非参数检验。我们以在执行需要正态分布数据的统计测试的许多情况下有用的不同数据转换方法结束本章。

在下一章中，我们将详细探讨假设检验，并讨论如何从测试结果中得出统计结论。我们还将探讨统计测试中可能出现的错误以及如何选择统计功效。

参考文献

• [1] 偏度 – www.itl.nist.gov/div898/handbook/eda/section3/eda35b.htm

• [2] 峰度 – [www.itl.nist.gov/div898/handbook/eda/section3/eda35b.htm#:~:text=峰度%20 是%20 一个%20 衡量,%20 会是%20 极端情况](https://www.itl.nist.gov/div898/handbook/eda/section3/eda35b.htm#:~:text=峰度%20 是%20 一个%20 衡量,%20 会是%20 极端情况).

• [3] 正态分布 – C.F. 高斯和最小二乘法，ŚLĄSKI PRZEGLĄD STATYSTYCZNY 西里西亚统计评论，第 12(18)期，O. Sheynin，1999 年 9 月

第三章：假设检验

在本章中，我们将开始讨论从数据中得出统计结论，结合来自第一章**，抽样和推广和来自第二章**，数据分布的分布。我们主要使用统计模型来回答数据中的问题。假设检验提供了一个正式的框架，通过不确定性的度量来回答感兴趣的问题。首先，我们将介绍假设检验的目标和结构。然后，我们将讨论假设检验可能出现的错误并定义预期的错误率。接着，我们将通过 z 检验的过程来讲解假设检验。最后，我们将讨论统计功效分析。

在本章中，我们将涵盖以下主要主题：

• 假设检验的目标

• 第一类和第二类错误

• z 检验的基础 – z 分数、z 统计量、临界值和 p 值

• 均值和比例的单样本和双样本 z 检验

• 选择错误率和功效分析

• 将功效分析应用于 z 检验

假设检验的目标

简单来说，假设检验的目标是判断我们所拥有的数据是否足以支持一个特定的假设。假设检验提供了一个正式的框架，基于我们的数据来检验假设，而不是试图通过视觉检查来做出决定。在本节中，我们将讨论假设检验的过程。在下一节中，z 检验的基础 – z 分数、z 统计量、临界值和 p 值，我们将通过详细分析一个例子来应用这个过程。

均值假设检验概述

为了理解假设检验过程，让我们从一个简单的例子开始。假设我们有一个工厂，工厂里有机器生产小部件，我们期望我们的机器以一定的速率（每小时 30 个小部件）生产小部件。我们首先构建两个假设，即零假设和备择假设。零假设和备择假设分别用以下符号表示：H 0 和 H a。为了创建零假设，我们将首先假设我们想要测试的是真实的。在我们的例子中，零假设将是机器的平均输出为每小时 30 个小部件。一旦我们确定了零假设，然后我们创建备择假设，它只是零假设的矛盾。在我们的例子中，备择假设是机器的平均输出不是每小时 30 个小部件。请注意，我们的假设没有指示任何方向性，也就是说，备择假设包含低于和高于预期值的值。这被称为双尾测试，意味着有两个备择假设。我们还有单尾测试。例如，如果我们说零假设是机器的平均输出超过每小时 30 个小部件，那么备择假设将是机器的平均输出低于每小时 30 个小部件。这一组假设将是一个单尾测试。

我们例子中的双尾零假设和备择假设可以用以下数学方式表述：

H 0 :  x = 30

H a :  x ≠ 30

一旦我们有了零假设和备择假设以及数据，我们将使用软件进行测试。现在我们先暂时搁置测试的实现细节（这些将在下一节中详细讨论）。统计测试有两种可能的结果：拒绝零假设或未能拒绝零假设。如果我们的机器的平均输出在零假设中声明的值上具有统计学上的差异，那么我们将拒绝零假设。这意味着，根据数据，零假设中声明的值不是平均值的合理值。然而，如果我们的机器的平均输出在零假设中列出的值上没有统计学上的差异，我们将未能拒绝零假设。这意味着，根据数据，零假设中声明的值是平均值的合理值。在运行零假设测试并得出结论后，我们将提供一个置信区间（下一节将讨论），并确定推理的范围。

推理范围

推理范围由讨论在第一章中的抽样设计决定。有两个问题需要考虑 – 总体是什么，以及样本是如何从总体中选择的？ 在这个例子中，让我们假设我们正在测试一个大工厂（可能是数百台机器）的机器的平均输出。那么，总体就是工厂中的机器。如果我们随机抽取机器样本，那么我们的结论可以推广到整个总体。

尽管我们的当前例子是现实的，但它相当简单。在其他情况下，可能会有额外的考虑。例如，工厂中的机器可能有不同的型号和不同的年龄，这可能会影响输出。在这种情况下，我们可以使用分层随机抽样并对每个层进行推断。

假设检验步骤

本节概述了从提出假设到得出结论的假设检验过程。随着我们继续本章的学习，请记住以下关键假设检验步骤：

• 陈述零假设和备择假设

• 执行统计测试

• 确定结论：拒绝或未能拒绝零假设

• 提供统计结论、置信区间和推理范围

这些步骤适用于任何统计测试，我们将继续遵循这一系列步骤进行假设检验。在下一节中，我们将讨论假设检验可能导致的错误类型。

第一类错误和第二类错误

虽然数据可以给我们一个关于分布特性的良好概念，但假设检验可能导致错误。错误可能发生，因为我们是从总体中抽取随机样本。虽然随机化使得样本包含抽样偏差的可能性降低，但并不能保证随机样本能够代表总体。假设检验可能导致两种可能的错误：

• 第一类错误：当零假设实际上为真时拒绝零假设

• 第二类错误：当零假设实际上为假时未能拒绝零假设

第一类错误

当假设检验导致拒绝零假设，但实际上零假设是正确的时，会发生第一类错误。例如，假设我们有一个数据分布，总体均值为 30。我们陈述的零假设为 H 0 :  _ x  = 30。我们为测试随机抽取样本，但样本中的随机值恰好位于分布的高侧。因此，测试结果建议我们应该拒绝零假设。在这种情况下，我们犯了一个第一类错误。这种类型的错误也称为假阳性。

当我们进行统计测试时，由于从目标总体中抽取的样本数据，我们总是有可能得出错误的结论。犯第一类错误的概率由 α指定。换句话说，α代表我们期望犯错误（预期错误率）的频率。这是一个我们可以为我们的测试选择的自由参数（α也称为显著性水平）。通常使用 0.05 作为α，但没有使用 0.05 的证据基础；在其他情况下，不同的值可能更合适。在本章的后面部分，我们将讨论选择第一类错误率。

第二类错误

我们可能犯的另一种错误称为第二类错误。在这种情况下，我们未能拒绝实际上为假的零假设。让我们考虑另一个例子。假设我们有一组数据分布，我们想要测试该分布的均值是否为 30。我们随机抽取样本进行测试，测试结果表明我们不应该拒绝零假设。然而，真正的总体均值是 35。在这种情况下，我们犯了一个第二类错误。这种错误也称为假阴性。

如前所述，统计测试可能导致错误的结论总是有可能的。因此，我们希望控制犯错误的概率。然而，与α不同，第二类错误率β不是一个我们可以简单选择的自由参数。为了理解犯第二类错误的可能性的大小，我们通常将进行功效分析，这将显示各种因素，如样本大小，将如何影响第二类错误率。在下一节中，我们将讨论选择错误率和功效分析。

我们可以用图 3.1中的表格来总结假设检验的可能结果。

		零假设是：
		真实
对零假设的决策	不拒绝	正确推断(1- α)
拒绝	第一类错误(α)	正确推断(1-β)

图 3.1 – 假设检验结果

在本节中，我们讨论了从统计测试中得出结论时可能出现的错误类型。在下一节中，我们将通过 z 检验的假设检验示例进行讲解，在本章的后面部分，我们将讨论如何选择错误率以及如何分析统计功效和相关因素。

z 检验的基础 – z 分数，z 统计量，临界值和 p 值

在本节中，我们将讨论一种称为 z 检验的假设检验类型。这是一种使用假设为正态分布的样本数据来确定与总体参数值相关的统计陈述是否应该被拒绝的统计程序。该测试可以在以下方面进行：

• 单样本（左尾 z 检验，右尾 z 检验或双尾 z 检验）

• 双样本（双样本 z 检验）

• 比例（单比例 z 检验或双比例 z 检验）

该测试假设标准差是已知的，并且样本量足够大。在实践中，样本量应大于 30。

在讨论不同类型的 z 检验之前，我们将讨论 z 分数和 z 统计量。

z 分数和 z 统计量

要衡量一个特定值与均值之间的距离，我们可以使用 z 分数或 z 统计量作为统计技术，结合均值和标准差来确定相对位置。

z 分数的计算公式如下：

z i =  x i −  _ x  _ σ

在这里，z i 是 x i 的 z 分数，_x 是样本均值，σ 是样本标准差。z 分数也称为 z 值、标准化值或标准分数。让我们考虑几个例子。标准差告诉我们样本与分布均值的距离。如果 z i = 1.8，那么这个点距离均值 1.8 个标准差。同样，如果 z i = − 1.5，那么这个点距离均值 1.5 个标准差。符号确定它是大于还是小于样本均值。z i 为 -1.5 表示小于均值，而 z i 为 1.8 表示大于均值。现在让我们通过一个例子来解释。在德克萨斯州达拉斯的一所高中（在美国），我们要求学生参加匿名智商测试，进行一些统计研究。从这所学校收集的数据呈正态分布，智商分数的总体均值 μ = 98，总体标准差 σ = 12。一名学生参加了智商测试，他的分数是 110。他的智商分数高于均值，但他想知道自己是否位于前 5%。首先，我们将使用 z 分数公式来计算它：

z student = 110 − 98 _ 12  = 12 _ 12 = 1.

学生可以检查 z 表（图 3*.2*），例如，从网站www.z-table.com获取 0.8413 的值。他位于前 1-0.8413 = 0.1587 或 15.87%的学校智商分数中。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_03_002.jpg

图 3.2 – z 表

在 Python 中，我们可以使用累积分布函数（CDF）来计算它：

import scipy
round(scipy.stats.norm.cdf(1),4)
# 0.8413

我们在 Python 中得到的值与 z 表检查中的值相同；在这个例子中，z 分数=1 时的值为 0.8413。

另一个例子是一个从智商调查中随机抽取的 10 个分数样本：

90, 78, 110, 110, 99, 115, 130, 100, 95, 93

要计算样本中每个智商分数的 z 分数，我们需要计算这个样本的均值和标准差，然后应用 z 分数公式。幸运的是，我们再次可以使用scipy库，如下所示：

import pandas as pd
import numpy as np
import scipy.stats as stats
IQ = np.array([90, 78,110, 110, 99, 115,130, 100, 95, 93])
z_score = stats.zscore(IQ)
# Create dataframe
data_zscore = {
  "IQ score": IQ,
  "z-score": z_score
}
IQ_zscore = pd.DataFrame(data_zscore)
IQ_zscore

我们创建了一个名为IQ的智商分数数组，并使用scipy.stats中的z-score来计算z_score。最后，我们创建了以下输出 DataFrame。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_03_003.jpg

图 3.3 – 输出 DataFrame

在讨论 z 统计量之前，我们将介绍抽样分布的概念。回顾上一个例子，作为一个经验法则，我们反复从高中智商分数池中选取 35 个智商分数的简单随机样本，所需次数取决于研究。然后，我们计算每个样本的平均分数，称为 x。因为我们有各种选出的样本，所以我们也有 x 的各种可能值。x 的期望值如下：

E(x) = μ

在这里，μ 表示总体均值。σ_x 表示 x 的标准分布。实际上，在许多抽样情况下，总体相对于较小的样本量来说相对较大。那么，x 的标准差可以表示如下：

σ_x = σ / √n

这里，σ 是总体标准差，n 是有限或无限总体中的样本量，总体较大，样本量相对较小。注意，σ_x 也被称为均值标准误，帮助我们确定样本均值与总体均值之间的距离。注意，E(x) = μ，与样本量无关。样本量和标准误成反比：当样本量增加时，标准误减小。由于假设 x 的抽样分布是正态分布的，样本分布由以下公式给出：

z_x = x̄ − μ / σ / x̄

在关于总体均值的假设检验中，我们使用测试统计量，其公式如下：

z = z_x = x − μ / σ / √n

在这里，x 是样本均值，μ 是总体均值，σ 是总体标准差，n 是样本量。

再次考虑智商测试的例子。智商分数数据的均值为 μ = 98，标准差 σ = 12。假设数据服从正态分布。设 x 为从智商数据中随机抽取的分数。x 在 95 和 104 之间的概率是多少？我们将计算当 x = 95 和 x = 104 时的 z 分数如下：

z_95 = 95 − 98 / 12 = −0.25,

z_100 = 104 − 98 / 12 = 0.5.

因此，所抽取的分数在 95 和 104 之间的概率是：

P(95 < x < 104) = P(−0.25 < z < 0.5) = 0.6915 − 0.4013 = 0.2902.

然后，这个概率大约是 29.02%，即从数据中随机抽取的智商分数在 95 和 104 之间的概率。我们还可以从 z 表中获取以下值。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_03_004.jpg

图 3.4 – z 表

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_03_005.jpg

图 3.5 – z 表

在 Python 中，我们可以按照以下方式实现代码：

#calculate z scores at x=95 and 104
zscore_95 = round((95-98)/12,2)
zscore_104 = round((104-98)/12,2)
#calculate cdf and probability
cdf_95 = stats.norm.cdf(zscore_95)
cdf_104 = stats.norm.cdf(zscore_104)
prob = abs(cdf_95-cdf_104)
#print the probability
print(f"The probability that the taken score between 95 and 104 is {round(prob*100,2)}%!")

提出的另一个问题是，随机抽取四个分数的平均分 x 在 95 和 104 之间的概率是多少？为了解决这个问题，我们使用 z 统计量的概念。我们假设 x 的均值也是 98，且 x 服从正态分布。那么，x 的标准误为：

σ_x = σ / √n = 12 / √4 = 6,

然后：

z =   _ x  − μ  _ x  _ σ  _ x   =   _ x  − 98 _ 6 .

很容易计算出 z 95 = − 0.5 和 z 104 = 1。通过使用 z 表，我们可以得到随机抽取六个分数的平均分数 x 在 95 和 104 之间的概率如下：

0.8413 – 0.3085 = 0.5328

或者，大约 53.28%。为了使用这个想法，我们可以在 Python 中实现以下代码：

import pandas as pd
import numpy as np
import scipy.stats as stats
import math
# standard error
n= 4
sigma = 12
se = sigma/math.sqrt(n)
#calculate z scores at x=95 and 104
zscore_95 = round((95-98)/se,2)
zscore_104 = round((104-98)/se,2)
#calculate cdf and probability
cdf_95 = stats.norm.cdf(zscore_95)
cdf_104 = stats.norm.cdf(zscore_104)
prob = abs(cdf_95-cdf_104)
#print the probability
print(f"The probability that the taken score between 95 and 104 is {round(prob*100,2)}%!")

注意，在上述代码中，我们还需要 math 库来计算平方根函数，math.sqrt()。

在下一节中，我们将讨论均值 z 测试。

均值 z 测试

在这部分，我们将考虑与总体均值或两个总体的均值相关的单样本和双样本 z 测试。

单样本 z 测试

从总体中选择的研究样本应该是正态分布的。总体标准差在实用目的上应该是已知的。

当总体不能假设为正态分布，但根据研究参与者的经验，样本大小需要足够大时，此测试仍然适用。为了进行假设检验，我们需要制定零假设和备择假设。以下图示了与左尾、右尾和双尾 z 测试相对应的三个零假设和备择假设。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_03_006.jpg

图 3.6 – 左尾、右尾和双尾假设检验

H 0 : μ ≥ μ 0 H 0 : μ ≤ μ 0 H 0 : μ = μ 0

H a : μ < μ 0 H a : μ > μ 0 H a : μ ≠ μ 0

接下来，我们需要指定显著性水平 α，即在零假设为真时拒绝零假设的概率。换句话说，这是我们在上一节中讨论的第一类错误的概率。然后，我们计算测试统计量的值。有两种方法用于假设检验：使用 p 值或临界值。

在 p 值方法中，我们使用测试统计量的值来计算一个概率，称为 p 值，它可能具有与从样本中得到的测试统计量一样极端或更极端的值。p 值越小，就越表明有证据反对零假设，换句话说，这是用来确定是否应该拒绝 H 0 的概率。拒绝规则（拒绝 H 0）是 p 值小于或等于研究中的指定显著性水平 α。为了根据测试统计量的值在 Python 中找到 p 值，我们使用以下语法：

 scipy.stats.norm.sf(abs(x))

在这里，x 是 z 分数。例如，我们想要找到一个左尾测试中与 z 分数 -2.67 相关的 p 值。Python 实现如下：

import scipy.stats
#find p-value
round(scipy.stats.norm.sf(abs(-2.67)),4)

输出将是 0.0038。在右尾测试的情况下，使用类似的 Python 代码。对于双尾测试，我们需要将值乘以 2：

#find p-value for two-tailed test
scipy.stats.norm.sf(abs(2.67))*2

下图说明了如何计算每种类型检验中的 p 值。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_03_007.jpg

图 3.7 – 假设检验中的 p 值

最后一步是解释统计结论。

另一方面，在临界值方法中，我们需要使用显著性水平来计算检验统计量的临界值。临界值是临界区域的边界，我们可以据此拒绝零假设。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_03_008.jpg

图 3.8 – 假设检验中的临界区域

在 Python 中计算临界值时，我们使用以下语法：

scipy.stats.norm.ppf(alpha)

这里，alpha 是要使用的显著性水平。以下是在 Python 中实现左尾、右尾和双尾检验的代码：

import scipy.stats
alpha = 0.05 # level of significance
#find Z critical value for left-tailed test
print(f" The critical value is {scipy.stats.norm.ppf(alpha)}")
#find Z critical value for left-tailed test
print(f" The critical value is {scipy.stats.norm.ppf(1-alpha)}")
##find Z critical value for two-tailed test
print(f" The critical values are {-scipy.stats.norm.ppf(1-alpha/2)} and {scipy.stats.norm.ppf(1-alpha/2)}")

这是前面代码的输出：

临界值 是 -1.6448536269514729

临界值 是 1.6448536269514722

临界值是 -1.959963984540054 和 1.959963984540054

对于左尾检验，显著性水平 α = 0.05，临界值约为 -1.64485。由于这是一个左尾检验，如果检验统计量小于或等于这个临界值，我们拒绝零假设。同样，对于右尾检验，如果检验统计量大于或等于 1.64485，我们拒绝零假设。对于双尾检验，如果检验统计量大于或等于 1.95996 或小于或等于 -1.95996，我们拒绝零假设。

确定是否拒绝零假设后，我们解释统计结论。

让我们再次讨论达拉斯高中的智商测试分数。智商分数数据平均值为 μ = 98，标准差为 σ = 12。一位研究人员想知道智商分数是否会受到某些智商训练的影响。他招募了 30 名学生，每天训练他们 2 小时回答智商问题，共训练 30 天，并在训练结束后记录他们的智商水平。训练阶段后 30 名学生的智商分数为 95, 110, 105, 120, 125, 110, 98, 90, 99, 100, 110, 112, 106, 92, 108, 97, 95, 99, 100, 100, 103, 125, 122, 110, 112, 102, 92, 97, 89, 和 102。他们的平均分数是 104.17。我们可以轻松地在 Python 中实现以下计算：

IQscores = [95,110, 105, 120, 125, 110, 98, 90, 99, 100,
            110, 112, 106, 92, 108, 97, 95, 99, 100, 100,
            103, 125, 122, 110, 112, 102, 92, 97, 89, 102]
IQmean = np.array(IQscores).mean()

我们定义零假设和备择假设：

H 0 : μ after = μ = 98

H a : μ after > μ = 98

我们选择显著性水平 α=0.05。之前，右尾检验的临界值为 1.64485。我们现在计算问题中的检验统计量：

z = x − μ σ/ √n = 104.17 − 98 12/ √30 = 2.8162

它在 Python 中的实现如下：

n=30 #number of students
sigma =12 #population standard deviation
IQmean = 104.17 # IQ mean of 30 students after the training
mu = 98 # population mean
z = (IQmean-mu)/(sigma/math.sqrt(n))

由于检验统计量值为 2.8162 > 1.64485，我们拒绝零假设。这意味着训练确实影响了这些学生的智商水平，并帮助他们提高了智商分数。

我们还可以使用来自statsmodels包的ztest()函数（Seabold, Skipper, 和 Josef Perktold, “statsmodels: Econometric and statistical modeling with python.” Proceedings of the 9th Python in Science Conference. 2010）来执行单样本或双样本 z 检验（我们将在下一部分讨论）。语法如下：

statsmodels.stats.weightstats.ztest(x1, x2=None,
    value=0, alternative='two-sided')

这里，我们看到以下内容：

• x1: 两个独立样本中的第一个

• x2: 两个独立样本中的第二个（如果执行双样本 z 检验）

• value: 在单样本情况下，value是零假设下x1的均值。在双样本情况下，value是零假设下x1的均值与x2的均值之差。检验统计量是x1_mean - x2_mean - value。

• alternative: 替代假设：

'two-sided': 双侧检验

'larger': 右侧检验

'smaller': 左侧检验

以下 Python 代码展示了我们如何执行单样本 z 检验：

from statsmodels.stats.weightstats import ztest as ztest
#IQ scores after training sections
IQscores = [95,110, 105, 120, 125, 110, 98, 90, 99, 100,
            110, 112, 106, 92, 108, 97, 95, 99, 100, 100,
            103, 125, 122, 110, 112, 102, 92, 97, 89, 102]
#perform one sample z-test
z_statistic, p_value = ztest(IQscores, value=98, alternative = 'larger')
print(f"The test statistic is {z_statistic} and the
    corresponding p-value is {p_value}.")

检验统计量为 3.3975，p 值=0.00034 < 0.05，其中 0.05 是显著性水平。因此，我们有足够的证据拒绝零假设。这意味着训练确实影响了这些学生的智商水平。

双样本 z 检验

我们考虑两个正态分布且独立的总体。让Ω表示两个总体均值μ1 和μ2 之间的假设差异。同样，正如单样本 z 检验的情况一样，我们有三种零假设和替代假设的形式：

H0: μ1 − μ2 ≥ Ω H0: μ1 − μ2 ≤ Ω H0: μ1 − μ2 = Ω

H0: μ1 − μ2 < Ω H0: μ1 − μ2 > Ω H0: μ1 − μ2 ≠ Ω

在许多问题中，Ω = 0。这意味着对于双尾检验的情况，零是零假设，换句话说，μ1 和μ2 是相等的。假设检验的检验统计量计算如下：

z = (̅x1 − ̅x2) − Ω / √(σ1²/n1 + σ2²/n2)

在这里，_x1 和 _x2 是样本均值，n1 和 n2 是从两个具有均值μ1 和μ2 的总体中随机抽取的样本大小。σ1 和σ2 是这两个总体的标准差。对于两个独立的简单随机样本，点估计量 _x1 − _x2 的标准误差如下：

σ̅x1−̅x2 = √(σ1²/n1 + σ2²/n2)

当样本量足够大时，_x1 − _x2 的分布可以被认为是正态分布。双样本 z 检验假设检验的逐步方法与单样本 z 检验类似。

让我们考虑一个例子。我们研究了达拉斯两所学校 A 和 B 的学生的智商分数，并想知道这两所学校的平均智商水平是否不同。记录了每所学校 30 名学生的简单随机样本。

A= [95,110, 105, 120, 125, 110, 98, 90, 99, 100,110, 112, 106, 92, 108, 97, 95, 99, 100, 100, 103, 125, 122, 110, 112, 102, 92, 97, 89, 102]

B = [98, 90, 100, 93, 91, 79, 90, 100, 121, 89, 101, 98, 75, 90, 95, 99, 100, 120, 121, 95,

96, 89, 115, 99, 95, 121, 122, 98, 97, 97]

研究人员收集了一个样本数据，观察到的样本比例 p = 0.84，假设的总体比例 p0 = 0.8，样本大小 n = 500。零假设和替代假设如下：

H0 : μ1 − μ2 = 0

Ha : _p ≠ p0.

我们选择显著性水平 α=0.05。在 Python 中，通过使用 statsmodels 包的 ztest() 函数，我们进行以下计算：

from statsmodels.stats.weightstats import ztest as ztest
#IQ score
A= [95,110, 105, 120, 125, 110, 98, 90, 99, 100,
    110, 112, 106, 92, 108, 97, 95, 99, 100, 100,
    103, 125, 122, 110, 112, 102, 92, 97, 89,102] #school A
B = [98, 90, 100, 93, 91, 79, 90, 100, 121, 89,
     101, 98, 75, 90, 95, 99, 100, 120, 121, 95,
     96, 89, 115, 99, 95, 121, 122, 98, 97, 97]  # school B
#perform two- sample z-test
z_statistic, p_value = ztest(A, B, value=0, alternative = 'two-sided')
print(f"The test statistic is {z_statistic} and the corresponding p-value is {p_value}.")

H0 : _p = p0,

比例的 z-test

H0 : _p ≥ p0 H0 : _p ≤ p0 H0 : _p = p0

单比例 z-test

在这里，n 是样本大小。让我们考虑一个例子。在休斯顿的一所社区学院，一位研究人员想知道学生是否支持一些等于 80% 的变化。他将在显著性水平 α = 0.05 下使用单比例 z-test。为了在 Python 中实现代码，他可以使用 statsmodel 库中的 proportions_ztest。语法如下：

nobs: 试验或观察的数量

单比例 z-test 用于比较样本比例 _p 与假设比例 p0 之间的差异。同样，正如在均值 z-test 中，我们对零假设和替代假设有三种形式 - 左尾、右尾和双尾测试：

测试统计量计算如下：

z = _p − p0___________√___________p0(1 − p0)_________n

在前面的代码中，我们选择了 alternative = 'two-sided'，这与研究的零假设和替代假设相关。Python 代码产生的测试统计量和 p-value 分别为 1.757 和 0.079。使用显著性水平 0.05，由于 p-value > 0.05，我们未能拒绝零假设。换句话说，我们没有足够的证据表明来自两所学校的学生的智商平均分数不同。

statsmodels.stats.proportion.proportions_ztest(count,
    nobs, value=None, alternative='two-sided')

• 我们也可以测试比例的差异。让我们看看如何进行比例的 z-test。

• 零假设和替代假设如下：

• value: 这是零假设的值，在单样本测试的情况下等于比例。在双样本测试的情况下，零假设是 prop[0] - prop[1] = value，其中 prop 是两个样本中的比例。如果没有提供，value = 0，零假设是 prop[0] = prop[1]。

  • 'two-sided': 双尾测试

  • 'larger': 右尾测试

  • 'smaller': 左尾测试

Ha : μ1 − μ2 ≠ 0

alternative: 替代假设：

Ha : _p < p0 Ha : _p > p0 Ha : _p ≠ p0

我们将在 Python 中实现一个单比例双尾 z-test 来计算测试统计量和 p-value：

#import proportions_ztest function
from statsmodels.stats.proportion import proportions_ztest
count = 0.8*500
nobs = 500
value = 0.84
#perform one proportion two-tailed z-test
z_statistic, p_value = proportions_ztest(count, nobs,
    value, alternative = 'two-sided')
print(f"The test statistic is {z_statistic} and the
    corresponding p-value is {p_value}.")

测试统计量为-2.236067977499786，相应的 p 值为 0.0253473186774685。由于 p 值 < 0.05（显著性水平），我们拒绝零假设。有足够的证据表明，支持变化的学生的比例与 0.8 不同。

双比例 z 检验

此检验用于检验两个总体比例之间的差异。还有三种形式的零假设和备择假设。测试统计量的计算如下：

z =  ‾ p 1 − ‾ p 2  ____________  √ _____________   _ p (1 −  _ p )( 1 _ n 1 +  1 _ n 2)

在这里，_ p 1 是从总体 1 的简单随机样本中的样本比例，_ p 2 是从总体 2 的简单随机样本中的样本比例，n 1 和 n 2 是样本大小，_ p  是总合并比例，计算如下：

_ p  =  n 1‾ p 1 + n 2‾ p 2 _ n 1 + n 2 .

如果 A 校支持变化的学生的比例与 B 校支持变化的学生的比例存在差异，我们将考虑一个类似的例子。在这里，n 1 = 100，n 2 = 100，_ p 1 = 0.8，_ p 2 = 0.7。你可以使用statsmodels库中的proportions_ztest函数（Seabold, Skipper, 和 Josef Perktold, “statsmodels: Econometric and statistical modeling with python.” Proceedings of the 9th Python in Science Conference. 2010）来进行假设检验。在这里，我们将直接使用测试统计量来计算 p 值。零假设和备择假设如下：

H 0 : _ p 1 = _ p 2 (两个总体比例相等)

H a : _ p 1 ≠ _ p 2 (两个总体比例不同)

接下来，我们指定双尾检验的显著性水平为α = 0.05。然后我们计算测试统计量和 p 值如下：

scipy.stats.norm.sf(abs(z))*2

测试统计量为 1.633，双尾检验的 p 值为 0.10。因为 p 值大于指定的显著性水平 0.05，所以我们未能拒绝零假设。有足够的证据表明，支持变化的 A 校和 B 校学生的比例是不同的。

最后，实现上述计算的 Python 代码如下：

import math
import scipy
p_1bar = 0.8
p_2bar = 0.7
n1 = 100.0
n2 = 100.0
p= (p_1bar*n1 + p_2bar*n2)/(n1+n2)   # the total pooled proportion
z = (p_1bar-p_2bar)/math.sqrt(p*(1-p)*(1/n1+1/n2))
pval = scipy.stats.norm.sf(abs(z))*2
print(f"The test statistic is {z} and the p-value for two tailed test is {pval}.")

在本节中，我们讨论了不同的重要统计概念，例如 z 分数、z 统计量、临界值、p 值以及均值和比例的 z 检验。我们将在下一节讨论选择错误率和功效分析。

选择错误率和功效分析

统计学概括了围绕形成从数据行为中得出可接受结论的近似。因此，统计错误是不可避免的。统计模型发现的重要性基本上由错误率决定。一种可以用来最小化建模错误的方法，尤其是在样本量较低的情况下，是功率测试。统计功率是正确拒绝零假设的概率，从而最小化 II 型错误。其中α（α）是 I 型错误，β（β）是 II 型错误，功率的公式是1 – β。为了回顾，I 型错误是错误拒绝零假设的概率。

正如所注，功率在样本量较小时更为重要，因为大数定律通常有助于最小化误差，前提是选择了适当的抽样方法。当比较的差异相对较小时，功率也很重要。一种功率分析可能特别有用的场景是抽样成本较高时，例如在人类研究中。功率分析可以用来根据所需的功率、I 型错误率和效应大小找到适当的最小样本量。这些参数之间的关系可以根据需要进行探索。效应大小是在假设检验中比较的数据之间的差异或相似性，例如平均值的标准化差异或相关性。在假设检验中，常见的显著性水平（I 型错误）通常在 0.01 到 0.1 之间，常见的功率水平在 0.8 到 0.9 之间，尽管这些指标是案例依赖的（https://www.ncbi.nlm.nih.gov/pmc/articles/PMC7745163/#:~:text=The%20ideal%20power%20of%20a,high%20as%200.8%20or%200.9）。

功率的一些显著特性如下：

• 必须有足够的功率来检测有意义的不同

• 功率随样本量增加

• 功率随效应大小增加

• 随着功率的增加，标准差（以及方差和标准误）减小

在以下图中，我们可以看到 I 型错误（α）是我们可能错误地假设零假设并得出没有统计显著差异的结论的地方，而实际上数据点不属于数据源 2（零假设），而是属于数据源 1。II 型错误（β）是我们可能错误地假设该区域的数据属于数据源 1（备择假设）的地方，而实际上它属于数据源 2。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_03_009.jpg

图 3.9 – 左尾双样本 z 检验中的误差可视化

如前图所示，我们看到了一个左尾 t 检验，将右侧的零分布与左侧的备择分布进行比较。在此图中，我们可以看到两个突出的概念：

• 随着α的减小，临界值向左滑动，更多地接近分布的异常值，β变大，反之亦然

• 统计功效是数据源 2 中β右侧的区域，即功效 = (1-β)

一类错误是基于预先选择的阈值确定的；研究者可能对 90%、95%或 99%的置信水平感到最舒适。然而，二类错误是基于参数化（标准差、效应量、样本量等）。因此，我们可以使用功效分析来确定样本量，反之亦然。各种假设检验的功效分析将在以下两章中实现。

双样本合并 z 检验的功效分析

让我们看看两个不同学科（学科 A 和学科 B）的教授薪资数据集。我们想知道根据我们拥有的数据，这两个学科之间是否存在统计上显著的差异。首先，我们需要进行功效分析，以了解我们是否有足够的样本来信任我们可能进行的 z 检验的结果。我们需要的功效分析组件包括效应量、一类错误率、期望的二类错误率、备择假设的方向（组 2 预计比组 1 大、小，或者可能比组 1 大或小），以及较大样本中观测值的比例相对于较小样本。学科 A 有 181 个薪资数据，学科 B 有 216 个。因此，216 将是我们的分子，对应于我们将考虑的组 1（学科 A 将是组 2）。

让我们假设我们不确定一个组是否会比另一个组大或小；我们将考虑这是一个双尾假设检验。这个 z 检验的效应量是两组之间的差异。我们将使用Cohen 的 d。为了计算它，我们需要将两组均值的差异除以合并标准差。在这里计算，合并标准差是组 1 样本数量乘以组 1 的方差，加上组 2 的相同值，所有这些除以两组的合并样本量。以下是合并标准差的方程：

√ _____________   n 1 σ 1 2 + n 2 σ 2 2  _ n 1 + n 2

对于效应量，我们需要将均值差除以效应量，如下所示：

|μ 1 − μ 2| ______________  √ ___________  n 1 σ 1 2 + n 2 σ 2 2  _ n 1 + n 2

注意，在效应量的方程中，我们取均值差的绝对值。这是因为，对于 Cohen 的 d，我们总是需要这个测试的正均值差。

如此计算，如果我们想要 80%的功效（II 型错误率是 20%）和 0.05 的 I 型错误率（显著性水平），我们需要在组 1 中 172 个样本和在组 2 中 145 个样本。然而，如果我们想要 90%的功效和 0.01 的显著性水平（99%置信度），我们需要在组 1 中 325 个样本和在组 2 中 274 个样本：

import statsmodels.api as sm
import math
df_prof = sm.datasets.get_rdataset("Salaries", "carData").data
df_prof_A = df_prof.loc[df_prof['discipline'] == 'A']
df_prof_B = df_prof.loc[df_prof['discipline'] == 'B']
def pooled_standard_deviation(dataset1, dataset2, column) -> float:
    pooledSD = math.sqrt(((len(dataset1) - 1)*(dataset1[column].std()**2)+(len(dataset2) - 1)*(dataset2[column].std()**2))/(len(dataset1) + len(dataset2) - 2))
    return pooledSD;
stdDeviation = pooled_standard_deviation(
    dataset1 = df_prof_A, dataset2=df_prof_B,
    column='salary')
from statsmodels.stats.power import NormalIndPower
effect = abs(df_prof_B['salary'].mean() –
    df_prof_A['salary'].mean() ) / stdDeviation
# The difference between two means divided by std if pooled 2-sample
alpha = 0.05
power = 0.8
ratio=1.19 # # of obs in sample 2 relative to sample 1
analysis = NormalIndPower()
result = analysis.solve_power(effect, power=power, nobs1=None, ratio=ratio, alpha=alpha, alternative='two-sided')
print('Sample Size Required in Sample 1: {:.3f}'.format(
    result*ratio)) # nobs1 is the sample size.
print('Sample Size Required in Sample 2: {:.3f}'.format(
    result)) # nobs2 is the sample size.

此代码的输出如下：

Sample Size Required in Sample 1: 171.620
Sample Size Required in Sample 2: 144.218
effect = abs(df_prof_B['salary'].mean() –
    df_prof_A['salary'].mean() ) / stdDeviation
alpha = 0.01
power = 0.9
ratio=1.19 # # of obs in sample 2 relative to sample 1
analysis = NormalIndPower()
result = analysis.solve_power(effect, power=power,
    nobs1=None, ratio=ratio, alpha=alpha,
    alternative='two-sided')
print('Sample Size Required in Sample 1: {:.3f}'.format(
    result*ratio)) # nobs1 is the sample size.
print('Sample Size Required in Sample 2: {:.3f}'.format(
    result)) # nobs2 is the sample size.

此代码的输出如下：

样本所需样本量 1: 325.346

样本所需样本量 2: 273.400

摘要

在本章中，我们介绍了假设检验的概念。我们从一个基本的假设检验轮廓开始，包含四个关键步骤：

• 陈述假设

• 执行测试

• 确定是否拒绝或未能拒绝零假设

• 在推断范围内得出统计结论

然后我们讨论了可能发生的潜在错误以及假阳性和假阴性，并定义了测试的预期错误率（alpha）和测试的效力（beta）。

我们还讨论了名为 z 检验的统计程序。这是一种使用假设为正态分布的样本数据进行假设检验的类型。在关于不同类型 z 检验的章节中，也介绍了 z 分数和 z 统计量，例如用于均值或比例的单样本或双样本 z 检验。

最后，我们讨论了功效分析的原理和动机，功效分析可以用来识别错误拒绝零假设的概率以及确定样本大小。我们还探讨了双总体合并 z 检验的分析参数。在这里，我们简要考察了效应量，这是我们在进行假设检验时寻找的影响值（治疗的效应）。我们将在下一章讨论功效分析，随着我们对假设检验不同应用的迭代。

在下一章中，我们将讨论更多的参数假设检验。虽然这些测试仍然需要分布假设，但假设将比 z 检验的假设不那么严格。

第四章：参数检验

在上一章中，我们介绍了假设检验的概念，并展示了 z 检验的几个应用。z 检验是称为参数检验的假设检验家族中的一类检验。参数检验是强大的假设检验，但参数检验的应用需要数据满足某些假设。虽然 z 检验是一个有用的检验，但它受到所需假设的限制。在本章中，我们将讨论几个更多的参数检验，这将扩展我们的参数工具集。更具体地说，我们将讨论 t 检验的各种应用，当存在多于两个数据子组时如何进行测试，以及皮尔逊相关系数的假设检验。我们将以参数检验的效力分析讨论结束本章。

在本章中，我们将涵盖以下主要主题：

• 参数检验的假设

• t 检验——一种参数假设检验

• 多于两个组的测试和方差分析（ANOVA）

• 皮尔逊 相关系数

• 效力分析示例

参数检验的假设

参数检验对人口数据做出假设，要求统计实践者在建模之前分析数据，尤其是在使用样本数据时，因为当真正的总体参数未知时，样本统计量被用作总体参数的估计。这些是参数假设检验的三个主要假设：

• 正态分布的人口数据

• 样本独立

• 等方差人口（在比较两个或更多组时）

在本章中，我们讨论了 z 检验、t 检验、ANOVA 和皮尔逊相关。这些测试用于连续数据。除了这些假设之外，皮尔逊相关要求数据包含成对样本。换句话说，每个被比较的组中必须有相等数量的样本，因为皮尔逊相关基于成对比较。

虽然这些假设是理想的，但有许多情况下这些假设无法得到保证。因此，了解这些假设具有一定的稳健性，这取决于测试是有用的。

正态分布的人口数据

因为在参数假设检验中，我们感兴趣的是对总体参数进行推断，例如均值或标准差，我们必须假设所选择的参数代表分布，并且可以假设数据有一个中心趋势。我们还必须假设统计量（从样本或抽样分布中取出的参数值）代表其相应的总体参数。因此，由于我们在参数假设检验中假设总体是正态分布的，样本也应该同样是正态分布的。否则，假设样本代表总体是不安全的。

参数假设检验高度依赖于均值，并假设它强烈地代表了数据的中心点（所有总体数据都围绕均值中心分布）。考虑两个分布的均值正在被比较以测试它们之间是否存在统计上显著的差异。如果分布是偏斜的，均值就不会是数据的中心点，因此不能很好地代表分布。由于这种情况，从比较均值得到的推断将不可靠。

对正态分布数据的稳健性

许多假设检验在用样本对总体进行估计时指定了自由度。自由度迫使模型假设所使用的分布中存在比实际存在的额外方差。尽管分析中的统计参数保持不变，但假设的额外方差使得集中趋势的度量更接近。换句话说，使用自由度迫使集中趋势的度量更中心地代表计算它们的分布。这是因为假设样本——虽然代表了它们的总体——但以误差范围代表它们的总体。因此，使用自由度的参数假设检验对违反正态分布数据的要求具有一定的稳健性。

在 图 4*.1 中显示的图表中，我们有一个略微偏斜的分布。一个应用了自由度，而另一个没有。我们可以看到，无论是否使用自由度，均值和中位数之间的距离都是相同的。然而，使用自由度的分布出现了更多的误差（更多的方差）：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_04_001.jpg

图 4.1 – 展示自由度的影响

当使用考虑均值的假设检验时，我们可以看到，均值，虽然不是中心（如中位数那样），当使用自由度时，相对于所有数据点，更接近分布的中心。由于参数假设检验使用均值作为中心点，这对于模型的有用性很重要，因为当使用自由度时，均值更能代表数据的中心点。这是对正态性有一定稳健性的主要原因。还有一些其他稳健性在于统计解释，例如在确定置信水平时；如果分布不是完美正态分布，例如，使用 90% 的置信水平而不是 99% 的置信水平可能更有益。

测试正态分布数据

确定一个分布是否为正态分布，从而可以用于参数假设测试，有多种方法。通常，对正态性的遵守程度由研究人员自行决定。本节中的方法在基于视觉检查以及应用统计显著性水平方面留下了一些讨论正态性的空间。

视觉检查

识别一个分布是否为正态分布的最佳测试方法是基于视觉检查。我们可以使用 分位数-分位数（QQ）图和直方图——以及其他测试方法——来对分布进行视觉检查。

在以下代码片段中，我们使用 scipy.stats 模块的 probplot 函数生成原始数据的图表以及 QQ 图：

import matplotlib.pyplot as plt
import scipy.stats as stats
import numpy as np
mu, sigma = 0, 1.1
normally_distributed = np.random.normal(mu, sigma, 1000)

在 图 4*.2 中，我们可以看到第一列是指数分布数据的直方图，在其下方是它的 QQ 图。由于点与接近 45 度红色线的近似程度非常远，该红色线代表纯正态分布，我们可以得出结论，数据不是正态分布的。通过视觉检查第二列中的数据，我们可以看到直方图显示了近似正态分布的数据集。这由其下方的 QQ 图得到证实，其中点主要近似于 45 度红色线。至于 QQ 图的尾部，这些数据点代表了偏斜度的密度。我们期望在正态分布的数据集中，大部分数据点将趋向于红色线的中心。在指数分布中，我们可以看到向左的重密度，红色线的左侧尾部，以及向右上方的线上的点稀疏分布。QQ 图可以从左到右阅读，与直方图中的分布相对应，其中最小值出现在 x 轴的左侧，最大值出现在右侧：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_04_002.jpg

图 4.2 – 使用 QQ 图和直方图评估正态性

QQ 图和直方图的视觉检查应该足以帮助研究人员得出正态性假设是否被违反的结论。然而，在可能不想进行视觉检查的情况下——例如在构建数据科学管道时——有其他方法可以提供关于正态性的具体测量。最常用的三种测试方法是 Kolmogorov-Smirnov、Anderson-Darling 和 Shapiro-Wilk 测试。

Kolmogorov-Smirnov 测试更关注数据的中心性。然而，如果数据在中心周围有较大的方差，则测试的效力会降低。Anderson-Darling 测试比中心更关注数据的尾部，如果数据有重尾和极端异常值，则更有可能识别出与正态性的不符。这两个测试在大样本量上表现良好，但样本量较小时效力并不强。我们考虑的第三个测试，Shapiro-Wilk，比 Kolmogorov-Smirnov 和 Anderson-Darling 测试更通用，因此对小样本量更稳健。基于这些特性，使用 Shapiro-Wilk 测试在自动化管道中可能更有用。或者，降低所应用测试的置信水平可能更好。

Kolmogorov-Smirnov

scipy.stats模块中的kstest函数，使用stats.norm.cdf（scipy的累积密度函数）执行这种单样本测试版本。双样本版本测试与指定的分布进行比较，以确定两个分布是否匹配。在双样本情况下，要测试的分布必须作为numpy数组提供，而不是下面代码片段中使用的stats.norm.cdf函数（图 4*.3*）。然而，这超出了测试正态性的范围，所以我们将不会探讨这一点。

Kolmogorov-Smirnov 将计算出的测试统计量与基于表格的临界值进行比较（kstest内部计算这个值）。与其他假设检验一样，如果测试统计量大于临界值，则可以拒绝给定分布是正态分布的零假设。如果 p 值足够低，也可以进行这种评估。测试统计量是给定分布中所有数据点与累积密度函数之间最大距离的绝对值。

Kolmogorov-Smirnov 特殊要求

Kolmogorov-Smirnov 测试要求数据围绕零中心并缩放到标准差为 1。所有数据都必须进行转换才能进行测试，但可以应用于预转换的分布；中心化和缩放后的分布不需要用于进一步的统计测试或分析。

在以下代码片段中，我们测试以确认正态分布数据集normally_distributed是否为正态分布。该数据集的平均值为 0，标准差为 1。输出确认数据是正态分布的。图 4*.3*中的图表显示了围绕均值为 0、标准差为 1 的正态分布，其右侧是相同分布的指数转换版本：

from scipy import stats
import numpy as np
mu, sigma = 0, 1
normally_distributed = np.random.normal(mu, sigma, 1000)

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_04_003.jpg

图 4.3 – 正态分布和指数分布数据

在这里，我们运行 Kolmogorov-Smirnov 测试：

stats.kstest(normally_distributed,
             stats.norm.cdf)

statsmodels的科尔莫哥洛夫-斯米尔诺夫检验为我们数据产生了以下结果：

KstestResult(statistic=0.0191570377833315, pvalue=0.849436919292824)

如果我们使用相同的数据，但将其指数变换为右偏态，相同的检验表明数据不再呈正态分布：

stats.kstest(np.exp(normally_distributed), stats.norm.cdf)

显著的 p 值确认了非正态性：

KstestResult(statistic=0.5375205782404135, pvalue=9.59979841227121e-271)

接下来，让我们取一个均值为 100，标准差为 2 的 1000 个样本的分布。我们需要将其中心化到均值为 0，单位方差（标准差为 1）。在下面的代码片段中，我们生成数据，然后执行缩放并将其保存到normally_distributed_scaled变量中：

mu, sigma = 100, 2
normally_distributed = np.random.normal(mu, sigma, 1000)
normally_distributed_scaled = (
    normally_distributed-normally_distributed.mean()) /
    normally_distributed.std()

现在数据已按要求中心化和缩放，我们使用科尔莫哥洛夫-斯米尔诺夫检验进行检查。正如预期的那样，数据被确认呈正态分布：

stats.kstest(normally_distributed_scaled, stats.norm.cdf)

这是输出：

KstestResult(statistic=0.02597307287070466, pvalue=0.5016041053535877)

安德森-达尔林

与科尔莫哥洛夫-斯米尔诺夫检验类似，scipy的anderson检验，我们可以测试其他分布，但默认参数指定正态分布，dist="norm"，假设零假设是给定的分布与正态分布在统计上是相同的。对于每个测试的分布，必须计算不同的临界值集。

安德森-达尔林与科尔莫哥洛夫-斯米尔诺夫比较

注意，虽然安德森-达尔林和科尔莫哥洛夫-斯米尔诺夫检验都使用累积密度频率分布来检验正态性，但安德森-达尔林检验与科尔莫哥洛夫-斯米尔诺夫检验不同，因为它对累积密度频率分布尾部的方差赋予比中间更大的权重。这是因为分布尾部的方差可以以比分布中间更小的增量来测量。因此，安德森-达尔林检验比科尔莫哥洛夫-斯米尔诺夫检验对尾部更敏感。与科尔莫哥洛夫-斯米尔诺夫检验一样，计算一个检验统计量，并将其与临界值进行比较。如果检验统计量大于临界值，则可以在指定的显著性水平上拒绝给定分布是正态分布的零假设。

这里，我们使用安德森-达尔林检验来检验一个均值为 19，标准差为 1.7 的随机正态概率分布。我们还测试了该数据的指数变换版本：

import matplotlib.pyplot as plt
import seaborn as sns
import numpy as np
mu, sigma = 19, 1.7
normally_distributed = np.random.normal(mu, sigma, 1000)
not_normally_distributed = np.exp(normally_distributed);

图 4*.4* 展示了数据的分布图：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_04_004.jpg

图 4.4 – 正态分布与重尾指数分布

在下面的代码和输出中，在图 4.5中，我们可以看到在所有显著性水平上分布都是正态分布的。回想一下，显著性水平是 p 值（即，显著性水平 = 15.0 表示 p 值为 0.15 或更小是显著的）：

from scipy import stats
import pandas as pd
import numpy as np
def anderson_test(data):
data = np.array(data)
    test_statistic, critical_values, significance_levels = stats.anderson(normally_distributed, dist='norm')
    df_anderson = pd.DataFrame({'Test Statistic':np.repeat(test_statistic, len(critical_values)), 'Critical Value':critical_values, 'Significance Level': significance_levels})
df_anderson.loc[df_anderson['Test Statistic'] >= df_anderson['Critical Value'], 'Normally Distributed'] = 'No'
    df_anderson.loc[df_anderson['Test Statistic'] <df_anderson['Critical Value'], 'Normally Distributed'] = 'Yes'
    return df_anderson;
mu, sigma = 19, 1.7
normally_distributed = np.random.normal(mu, sigma, 1000)
anderson_test(normally_distributed)

在这里，通过 numpy 的 random.normal 函数生成数据，使用 Anderson-Darling 方法进行测试，并确认其为正态分布：

检验统计量	临界值	显著性水平	正态分布
0.191482344	0.574	15	是
0.191482344	0.653	10	是
0.191482344	0.784	5	是
0.191482344	0.914	2.5	是
0.191482344	1.088	1	是

图 4.5 – 正态分布数据的 Anderson-Darling 结果

在这里，我们测试了正态分布数据的指数变换，以检查其正态性。数据呈指数分布，应在所有显著性水平上被拒绝。然而，我们在图 4.6中看到，它在 0.01 的显著性水平（99% 置信度）上未能被拒绝。因此，根据具体的使用情况，检查所有显著性水平、使用不同的测试或基于多个测试的结果做出决定可能是谨慎的：

not_normally_distributed = np.exp(normally_distributed)
anderson_test(not_normally_distributed)

我们的非正态分布数据的 Anderson-Darling 测试结果如下所示在图 4.6：

检验统计量	临界值	显著性水平	正态分布
0.96277351	0.574	15	否
0.96277351	0.653	10	否
0.96277351	0.784	5	否
0.96277351	0.914	2.5	否
0.96277351	1.088	1	是

图 4.6 – 非正态分布数据的 Anderson-Darling 结果

Shapiro-Wilk

scipy.stats shapiro 模块，因此输入数据在测试之前不需要更改。在 scipy 中，此测试的显著性水平为 0.05。

Shapiro-Wilk 与 Kolmogorov-Smirnov 和 Anderson-Darling 的比较

与 Kolmogorov-Smirnov 和 Anderson-Darling 相比，Shapiro-Wilk 对于测试大约小于 50 的较小样本量是理想的。然而，一个缺点是，由于 Shapiro-Wilk 使用重复抽样和测试来通过蒙特卡洛模拟应用计算出的检验统计量，大数定律的风险在于，随着样本量的增加，存在固有的风险增加，遇到第二类错误（功率损失）并未能拒绝零假设，其中零假设表示给定的分布是正态分布的。

使用与 Anderson-Darling 测试相同的分布，我们使用 Shapiro-Wilk 进行测试。我们可以看到，具有均值为 19 和标准差为 1.7 的随机正态分布，Shapiro-Wilk 测试已确认，p 值为 0.99，不应拒绝输入分布是正态分布的零假设：

mu, sigma = 19, 1.7
normally_distributed = np.random.normal(mu, sigma, 1000)
stats.shapiro(normally_distributed)

这是输出结果：

ShapiroResult(statistic=0.9993802905082703, pvalue=0.9900037050247192)

当使用正态分布数据的指数变换版本进行测试时，我们发现显著的 p 值（p = 0.0），这表明我们有足够的证据拒绝零假设，并得出结论分布不是正态分布：

not_normally_distributed = np.exp(normally_distributed)
stats.shapiro(not_normally_distributed)

这是输出：

ShapiroResult(statistic=0.37320804595947266, pvalue=0.0)

独立样本

在参数假设检验中，样本的独立性是另一个重要假设。非独立抽样可能导致两种影响。一种影响发生在进行子群抽样时。这里的问题是，总体中一个子群的响应可能与另一个子群的响应不同，甚至可能更接近于不同总体的响应。然而，当抽样代表总体时，这种子群差异可能并不非常具有代表性。

非独立抽样的另一个影响是，当样本在时间上足够接近，以至于一个样本的发生会阻止或排除另一个样本的发生时。这被称为序列（或自）相关性。

参数检验通常对违反此要求不具鲁棒性，因为它对测试结果的可解释性有直接、分类的影响。关于子群抽样，可以通过诸如在第一章**，抽样与泛化中概述的良好结构化抽样方法来防止这种情况。然而，关于序列效应，我们可以在数据中测试自回归相关性（也称为序列相关性）。

杜宾-沃森

评估抽样独立性缺失的最常见测试之一是一阶自回归测试，称为杜宾-沃森测试。自回归意味着使用先前数据点来预测当前数据点。一阶意味着最后一个抽样数据点（滞后一）是与序列中最近抽样数据点（滞后零）最显著相关的点。在一阶自相关中，每个数据点的相关性最强的是前一个数据点。杜宾-沃森测试并不测试任何值是否与它之前的值相关，而是测试每个值与它之前值之间是否存在足够强的关系，以至于可以得出存在显著自相关性的结论。从这个意义上说，对非独立抽样有一定的鲁棒性，即一次或两次意外事件可能不会完全使假设检验无效，但此类违规的持续发生将会。

杜宾-沃森值为 2 表示没有显著的自相关性，介于 0 到 2 之间的值表示正（直接）自相关性，而介于 2 到 4 之间的值表示负（反向）自相关性。

在以下示例中，我们有两个分布，每个分布都有 1,000 个样本。左边的分布是一个正弦分布，表现出强烈的自回归相关性，而右边的分布是一组随机生成数据，显示出白噪声方差（围绕均值为 0 的随机点）。使用statsmodels.stats模块中的durbin_watson()函数，我们能够确认正弦模式中的直接、正的一阶自相关性（非常小的 Durbin-Watson 值）和与随机噪声的 Durbin-Watson 统计量为 2.1，表明没有自相关性。因此，在图 4.7 中，左边的图不是由独立样本组成，而右边的图是：

from statsmodels.stats.stattools import durbin_watson
import matplotlib.pyplot as plt
import numpy as np
mu, sigma = 0, 1.1
independent_samples = np.random.normal(mu, sigma, 1000)
correlated_samples = np.linspace(-np.pi, np.pi, num=1000)
fig, ax = plt.subplots(1,2, figsize=(10,5))
ax[0].plot(correlated_samples, np.sin(correlated_samples))
ax[0].set_title('Durbin Watson = {}'.format(
    durbin_watson(correlated_samples)))
ax[1].plot(independent_samples)
ax[1].set_title('Durbin Watson = {}'.format(
    durbin_watson(independent_samples)))

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_04_007.jpg

图 4.7 – 串联相关和正态分布的序列数据

等样本方差

与正态分布数据的假设类似，等样本方差的假设——也称为方差齐性——是关于比较的分布的物理属性形状。假设等样本方差有助于提高参数测试的功效。这是因为当均值被识别为不同时，我们有信心；我们还知道潜在分布重叠的程度。当测试对整个分布的位置有直观认识——实际上，对效应大小的真正了解——功效会增加。相反，随着方差的发散，功效会降低。

对等样本方差的鲁棒性

虽然等样本方差在参数测试中很有用，但存在对这些测试的修改，有助于使结果对偏离等方差的偏差具有鲁棒性。这些测试的一个突出修改版本是使用Welch-Satterthwaite自由度调整。因为当每个组有不同的方差时，将相同的自由度应用于每个组会导致数据的误表示，Welch-Satterthwaite 调整在分配假设等方差的参数测试的自由度时考虑了方差差异。使用 Welch-Satterthwaite 调整的两个常见测试是 Welch 的 t 检验和 Welch 的方差分析测试。当用于小样本时，这些测试可能不可靠，但当用于足够大的样本量以具有足够功效时，结果应与它们的非 Welch 对应物大致相同。

测试等方差

在测试分布之间的方差相等时，我们有两种突出的测试：Levene 方差齐性检验和Fisher 的 F 检验。

Levene 方差齐性检验

Levene 方差齐性检验在测试两个或更多组的方差齐性时很有用。在下面的代码片段中 图 4*.8*，我们使用三个分布进行测试，每个分布的样本量为 100，均值为 0，标准差为 0.9、1.1 和 2。图 4*.8* 是使用上述 图 4*.8* 代码输出的数据生成的三个分布的图。

from scipy.stats import levene
np.random.seed(26)
mu1, sigma1, mu2, sigma2, mu3, sigma3 = 0,0.9,0,1.1,0,2
distro1, distro2, distro3 = pd.DataFrame(), pd.DataFrame(),
    pd.DataFrame()
distro1['x'] = np.random.normal(mu1, sigma1, 100)
distro2['x'] = np.random.normal(mu2, sigma2, 100)
distro3['x'] = np.random.normal(mu3, sigma3, 100)

我们可以看到它们不同的标准差如何影响它们的范围。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_04_008.jpg

图 4.8 – 多重方差齐性检验的分布

我们可以看到，由于这个结果是一个具有统计学意义的 p 值，表明分布不具有同质性，因此测试对非同质性方差的违反很敏感：

f_statistic, p_value = levene(distro1['x'], distro2['x'], distro3['x'])
if p_value <= 0.05:
    print('The distributions do not have homogenous variance. P-value = %.4f, F-statistic = %.4f'%(p_value, f_statistic))
else:
    print('The distributions have homogenous variance.P-value = %.4f, F-statistic = %.4f'%(p_value, f_statistic))

这是输出结果：

分布不具有同质性。P 值 = 0.0000

Fisher 的 F 检验

当同时测试两组的方差齐性时，Fisher 的 F 检验很有用。这个检验将一个检验统计量与临界值进行比较，以确定方差是否在统计学上相同。计算出的 F 统计量是第一组的方差除以第二组的方差。第一组总是具有较大方差的组。使用前面的数据，让我们比较分布 1 和分布 3。分布 3 的方差较大，为 2，因此该组的方差将是计算 F 统计量时的分子。由于每个组都有 100 个样本量，它们的自由度在查表时各为 99。然而，由于我们将使用 scipy Python 包来计算检验，这里不需要查表，因为 scipy 使用 f.cdf() 函数为我们完成这项工作。与 Levene 检验的结果一致，F 检验表明分布 1 和分布 3 不具有同质性：

from scipy.stats import f
def f_test(inputA, inputB):
    group1 = np.array(inputA)
    group2 = np.array(inputB)
    if np.var(group1) > np.var(group2):
        f_statistic = np.var(group1) / np.var(group2)
        numeratorDegreesOfFreedom = group1.shape[0] - 1
        denominatorDegreesOfFreedom = group2.shape[0] - 1
    else:
        f_statistic = np.var(group2)/np.var(group1)
        numeratorDegreesOfFreedom = group2.shape[0] - 1
        denominatorDegreesOfFreedom = group1.shape[0] - 1
    p_value = 1 - f.cdf(f_statistic,numeratorDegreesOfFreedom, denominatorDegreesOfFreedom)
    if p_value <= 0.05:
print('The distributions do not have homogenous variance. P-value = %.4f, F-statistic = %.4f'%(p_value, f_statistic))
    else:
        print('The distributions have homogenous variance. P-value = %.4f, F-statistic = %.4f'%(p_value, f_statistic))
f_test(distro3['x'], distro1['x'])

这个 F 检验输出如下：

分布不具有同质性。P 值 = 0.0000，F 统计量 = 102622.9745

t 检验 – 参数假设检验

在上一章中，当总体标准差已知时，应用了均值 z 检验。然而，在现实世界中，获得总体标准差并不容易（或几乎不可能）。在本节中，我们将讨论另一种称为 t 检验的假设检验，当总体标准差未知时使用。通过取代表性总体样本数据的均值和标准差来估计总体均值和标准差。

从广义上讲，均值 t 检验的方法与均值 z 检验的方法非常相似，但检验统计量和 p 值的计算方法与 z 检验不同。检验统计量通过以下公式计算：

t = _x − μ _ s/ √ _ n

在这里，x，μ，s 和 n 分别是样本均值、总体均值、样本标准差和样本大小，当样本数据 x 正态分布时，它具有 t**-分布**。以下代码展示了标准正态分布（蓝色曲线）和两个样本大小（3 和 16 个样本）的 t-分布（绿色和红色曲线）：

# libraries
import numpy as np
import scipy.stats as stats
# creating normal distribution
x =np.linspace(-5, 5, 1000) #create 1000 point from -5 to 5
y = stats.norm.pdf(x) # create probability density for each point x  - normal distribution
# creating Student t distributions for 2 sample sizes n =3 and n =15
degree_freedom1 = 2
t_dis1 = stats.t.pdf(x, degree_freedom1)
degree_freedom2 = 15
t_dis2 = stats.t.pdf(x, degree_freedom2)

以下可视化是针对考虑的这 3 个分布。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_04_009.jpg

图 4.9 – 正态分布和 t-分布

观察到三条曲线具有相似的对称性和形状，但样本大小较小时，变异性更大（或者说，尾部更重）。历史上，当样本大小大于 30 时，研究人员认为样本标准差可以代表总体，即当 n > 30 时，红色曲线近似于蓝色曲线。如果样本分布与标准正态分布重叠，也常用 z-检验。这种做法背后有一些理由，因为以前，临界值表存储到样本大小为 50，但现在，随着计算能力和互联网的强大，任何样本大小的 t 值都可以轻松获得。

均值 t-检验

在这部分将考虑与总体均值或两个总体均值相关的单样本和双样本 t-检验，其中总体方差或总体标准差是未知的。

要进行 t-检验，需要满足以下假设：

正态性：样本是正态分布的

独立性：观察值是从总体中随机选取以形成样本，换句话说，它们是独立的

让我们在下一节考虑单样本 t-检验。

单样本 t-检验

与单样本 z-检验类似，为了进行假设检验，需要考虑零假设和备择假设。以下列出了对应于左尾、右尾和双尾 t-检验的三个零假设和备择假设：

H 0 : μ ≥ μ 0 H 0 : μ ≤ μ 0 H 0 : μ = μ 0

H a : μ < μ 0 H a : μ > μ 0 H a : μ ≠ μ 0

接下来，需要根据研究目的指定显著性水平 α。有两种方法：p 值方法和临界值方法。在 p 值方法中，拒绝规则（拒绝 H 0——零假设）是当 p 值小于或等于所选的指定显著性水平时。在临界值方法中，拒绝规则是当检验统计量小于或等于左尾 t-检验的临界值 - t α，对于右尾 t-检验，检验统计量大于或等于 t α，对于双尾检验，检验统计量小于或等于 - t α/2 或大于或等于 t α/2。最后一步是对假设检验的统计结论进行解释。

要根据学生 t 分布的值找到 p 值，我们可以使用以下语法：

scipy.stats.t.sf(abs(x), df)

在这里，x 是检验统计量，df 是自由度（df = n-1，其中 n 是样本大小）在公式中。

例如，为了找到一个左尾检验中自由度为 14 的 t 值为 1.9 的 p 值，这将是 Python 的实现：

import scipy.stats
round(scipy.stats.t.sf(abs(1.9), df=14),4)

输出将是 0.0391。如果显著性水平 α = 0.05，那么我们拒绝零假设，因为 p 值小于 α。对于右尾 t 检验，与左尾 t 检验类似的 Python 代码被实现以找到 p 值。对于双尾检验，我们需要将值乘以 2，如下所示：

scipy.stats.t.sf(abs(t), df)*2

在这里，t 是检验统计量，df 是自由度（df = n-1，其中 n 是样本大小）。

要在 Python 中计算临界值，我们使用以下语法：

scipy.stats.t.ppf(q, df)

在这里，q 是显著性水平，df 是公式中要使用的自由度。以下是 Python 中左尾、右尾和双尾检验的代码实现：

import scipy.stats as stats
alpha = 0.05 # level of significance
df= 15 # degree of freedom
#find t critical value for left-tailed test
print(f" The critical value is {stats.t.ppf(q= alpha, df =df)}")
#find t critical value for right-tailed test
print(f" The critical value is {stats.t.ppf(q= 1-alpha, df =df)}")
##find t critical value for two-tailed test
print(f" The critical values are {-stats.t.ppf(q= 1-alpha/2, df =df)} and {stats.t.ppf(q= 1-alpha/2, df =df)}")

这是前面代码的输出：

• 临界值为 -1.7530503556925552

• 临界值为 1.7530503556925547

• 临界值为 -2.131449545559323 和 2.131449545559323

在显著性水平 α = 0.05 下，对于左尾检验，临界值约为 -1.753。由于这是一个左尾检验，如果检验统计量小于或等于这个临界值，我们拒绝零假设。同样，对于右尾检验，如果检验统计量大于或等于 1.753，我们拒绝零假设。对于双尾检验，如果检验统计量大于或等于 2.1314 或小于或等于 -2.1314，我们拒绝零假设。最后，我们解释假设检验的统计结论。

让我们从一所高中随机选择 30 名学生并评估他们的智商。我们想测试这个高中学生智商分布的平均值是否高于 100。这意味着我们将执行右尾 t 检验。30 名学生的智商分数如下：

IQscores = [113, 107, 106, 115, 103, 103, 107, 102, 108, 107, 104, 104, 99, 102, 102, 105, 109, 97, 109, 103, 103, 100, 97, 107,116, 117, 105, 107, 104, 107]

在进行假设检验之前，我们将检查正态性和独立性假设。如果样本是从该学校高中学生的总体中随机选择的，则满足独立性假设。对于正态性，我们将检查智商分数的直方图和 QQ 图：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_04_010.jpg

图 4.10 – 视觉评估学生智商分数的正态性

从直方图和 QQ 图来看，几乎没有证据表明该高中学生智商分布的总体不是正态分布。由于假设分布是正态的，我们将继续进行 t 检验。

首先，我们定义零假设和备择假设：

H 0 : μ ≤ 100

H a : μ > 100

我们选择显著性水平α=0.05。你可以通过手动使用其数学公式或通过实现 Python 来计算检验统计量。对于临界值和 p 值，之前部分展示的实现代码如下。在这里，我们将使用scipy库中的另一个函数来找到检验统计量和 p 值：

scipy.stats.ttest_1samp(data, popmean, alternative='greater')

这里适用以下规则：

• data: 样本中的观测值

• popmean: 零假设中的期望值

• alternative: 对于双尾 t 检验为'two-sided'，对于左尾 t 检验为'less'，对于右尾 t 检验为'greater'

Python 代码实现如下：

import scipy.stats as stats
#perform one sample t-test
t_statistic, p_value = stats.ttest_1samp(IQscores, popmean =100, axis=0,  alternative='greater')
print(f"The test statistic is {t_statistic} and the corresponding p-value is {p_value}.")

这是输出：

检验统计量为 6.159178830896832，相应的 p 值为 5.15076734562176e-07。

因为 p 值<0.05，其中 0.05 是显著性水平，我们有足够的证据拒绝零假设，并得出结论，该学校学生的真实平均 IQ 分数高于 100。

此外，以 95%的置信度，平均 IQ 分数介于 104.08 和 107.12 之间。我们可以在 Python 中如下进行置信区间的计算：

IQmean = np.array(IQscores).mean() # sample mean
IQsd = np.array(IQscores).std() # sample standard deviation
sample_size = len(np.array(IQscores)) # sample size
df = sample_size-1 # degree of freedom
alpha = 0.05 # level of significance
t_crit = stats.t.ppf(q=1-alpha, df =df) # critical
confidence_interval = (IQmean-IQsd*t_crit/np.sqrt(sample_size), IQmean+IQsd*t_crit/np.sqrt(sample_size))

在 Python 中使用单尾 t 检验进行假设检验的步骤与右尾和双尾 t 检验的步骤相似。

双样本 t 检验 – 池化 t 检验

与第三章**假设检验（均值的双样本 z 检验）中所述类似，均值的双样本 t 检验对于零假设和备择假设有三种形式。在进行测试之前，需要满足一些假设，如下所示：

• 正态性：两个样本是从它们的正态分布总体中抽取的

• 独立性：一个样本的观测值相互独立

• 方差齐性：两个总体假设具有相似的标准差

对于正态性，我们使用视觉直方图和两个样本的 QQ 图进行比较。我们假设独立性得到满足。为了检查两个样本的方差是否相等，我们可以通过观察它们的直方图进行可视化，如果可视化结果不明确，还可以使用 F 检验来获得额外的证据。这是一个假设检验，用于检查两个样本方差是否相等。

让我们看看两所高中 A 和 B 之间的 IQ 分数。以下是从每所学校随机选取的 30 名学生的分数：

IQscoresA=[113, 107, 106, 115, 103, 103, 107, 102,108, 107,
            104, 104, 99, 102, 102, 105, 109, 97, 109, 103,
            103, 100, 97, 107, 116, 117, 105, 107, 104, 107]
IQscoresB = [102, 108, 110, 101, 98, 98, 97, 102, 102, 103,
             100, 99, 97, 97, 94, 100, 104, 98, 92, 104,
            98, 95, 92, 111, 102, 112, 100, 103, 103, 100]

图 4.11 中显示的直方图和 QQ 图是由上述 IQ 数据生成的。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_04_011.jpg

图 4.11 – 评估两所学校 IQ 分数的正态性

我们可以看到正态性假设得到满足。我们还可以通过观察直方图假设方差齐性假设得到支持。如果需要，以下是一个额外的 F 检验来检查方差齐性假设：

# F-test
import numpy as np
import scipy.stats as stats
IQscoresA = np.array(IQscoresA)
IQscoresB = np.array(IQscoresB)
f = np.var(IQscoresA, ddof=1)/np.var(IQscoresB, ddof=1) # F statistic
dfA = IQscoresA.size-1 #degrees of freedom A
dfB = IQscoresB.size-1 #degrees of freedom B
p = 1-stats.f.cdf(f, dfA, dfB) #p-value

前述代码的输出告诉我们 F 检验统计量是 0.9963，相应的 p 值是 0.50394 > 0.05（0.05 是显著性水平），因此我们未能拒绝零假设。这意味着有足够的证据表明这两个样本的标准差是相等的。

我们现在定义零假设和备择假设：

H 0 : μ A = μ B,

H a : μ A ≠ μ B.

我们选择显著性水平 α=0.05。我们使用 statsmodels.stats.weightstats.ttest_ind 函数进行 t 检验。文档可以在这里找到：www.statsmodels.org/dev/generated/statsmodels.stats.weightstats.ttest_ind.xhtml。

我们可以使用此函数通过 alternative='two-sided'，'larger' 或 'smaller' 来执行备择假设的三种形式。在合并方差 t 检验中，当满足方差相等的假设时，检验统计量的计算如下：

t =  (‾ x 1 − ‾ x 2) − (μ 1 − μ 2)  _____________  s p √ _  1 _ n 1 +  1 _ n 2

在这里，_ x 1, _ x 2, μ 1, μ 2, n 1, 和 n 2 分别是两个样本 1 和 2 的样本均值、总体均值和样本大小，这里给出了合并标准差：

s p = √ ________________   (n 1 − 1) s 1 2 + (n 2 − 1) s 2 2  ________________  n 1 + n 2 − 2

自由度在这里显示：

df = n 1 + n 2 − 2.

让我们回到例子：

from statsmodels.stats.weightstats import ttest_ind as ttest
t_statistic, p_value, degree_freedom = ttest(IQscoresA,
    IQscoresB, alternative='two-sided', usevar='pooled')

输出返回的检验统计量是 3.78，p 值是 0.00037，t 检验中使用的自由度是 58（每个样本大小有 30 个观测值，因此自由度计算为 30 + 30 - 2 = 58）。

因为 p 值 <0.05，我们拒绝零假设。有足够的证据表明，高中 A 和 B 的学生平均智商分数之间存在差异。为了进行置信水平测试，你可以适应单样本 t 检验最后一部分中的 Python 代码。

回想一下，在某些情况下，如果直方图和 QQ 图显示一些偏斜的证据，我们可以考虑在假设检验中测试中位数而不是均值。如第二章**，数据分布所示，我们可以执行数据转换（例如，对数转换）以获得正态性假设。转换后，log(data) 的中位数等于 log(data) 的均值。这意味着测试是在转换数据的均值上进行的。

双样本 t 检验 – Welch 的 t 检验

这是一个实际的双样本 t 检验，当数据呈正态分布但总体标准差未知且不相等时。我们对于正态性的假设与合并 t 检验相同，但在进行 Welch 的 t 检验时，我们可以放宽方差相等的假设。让我们考虑以下例子，其中我们有两个样本数据集：

sample1 = np.array([2,3,4,2,3,4,2,3,5,8,7,10])
sample2 = np.array([30,26,32,34,28,29,31,35,36,33,32,27])

我们假设独立性得到满足，但我们将检查这两个样本的正态性和等方差假设，如下所示：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_04_012.jpg

图 4.12 – 检查 Welch 的 t 检验的等方差

通过查看图 4.12中直方图的x轴刻度，可以强烈地看出标准差不等的证据。让我们假设正态性假设得到满足。在这种情况下，进行双样本合并 t 检验不是一个好主意，但 Welch 的 t 检验就足够了。这里给出了零假设和备择假设：

H 0 :  μ 1 = μ 2

H a :  μ 1 ≠ μ 2

我们指定显著性水平为 0.05。为了计算检验统计量和 p 值，我们实现以下代码：

import scipy.stats as stats
t_statistic, p_value = stats.ttest_ind(sample1, sample2,
    equal_var = False)

检验统计量为-22.47，p 值<0.05（显著性水平）。我们拒绝零假设。有强有力的证据表明样本数据 1 的均值与样本数据 2 的均值不同。

配对 t 检验

配对 t 检验也称为匹配对或相关 t 检验，在研究样本中的每个元素被测试两次（前测和后测或重复测量）且研究人员认为存在某些相似性（如家庭）时使用。这里列出了假设：

• 差异呈正态分布

• 差异在观察之间是独立的，但从一个测试到另一个测试是相关的

配对 t 检验在许多研究中被使用，尤其是在与前后治疗相关的医学推理测试中。让我们回顾一下智商测试分数——一位研究人员招募了一些学生，以查看训练前后是否有分数差异，如下表所示：

学生	训练前分数	训练后分数	差异
A	95	95	0
B	98	110	12
C	90	97	7
D	115	112	-3
E	112	117	5

图 4.13 – 训练前和训练后分数

在这种情况下，我们不应使用独立的双样本 t 检验。应该检验差异的均值。我们可以通过直方图和 QQ 图来检查关于正态分布的假设，如下所示：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_04_014.jpg

图 4.14 – 检查配对 t 检验的正态性

通过查看 QQ 图，比直方图更能明显地看出数据呈正态分布的证据。

假设差异是独立的。这里给出了零假设和备择假设：

H 0 : μ pos − μ pre = 0

H a : μ pos − μ pre > 0

d i 表示每个学生的训练前分数和训练后分数之间的差异。零假设和备择假设可以重写如下：

H 0 : μ d = 0

H a : μ d > 0

然后，计算检验统计量如下：

t = _ d  − μ d _  s d _ √ _ n

在这里，d 是差异的样本均值，s_d 是差异的样本标准差。换句话说，配对 t 检验被简化为单样本 t 检验。然而，我们可以直接在 scipy 中使用以下函数：

stats.ttest_rel(data_pos, data_pre, alternative = {'two-sided', '``less', 'greater'})

备择假设对应于左尾、右尾或双尾测试。以下是 IQ 测试分数研究示例的 Python 实现：

from scipy import stats
IQ_pre = [95, 98, 90, 115, 112]
IQ_pos = [95, 110, 97, 112, 117]
t_statistic, p_value = stats.ttest_rel(IQ_pos, IQ_pre, alternative = 'greater')

测试统计量是 1.594，p 值是 0.093。因此，鉴于 p 值<0.05 和显著性水平α=0.05，我们拒绝零假设。有足够的证据表明，培训对智商分数有显著影响。

多组测试和方差分析

在前一章和前几节中，我们介绍了两组之间的测试。在本节中，我们将介绍两种测试组间差异的方法，如下：

• 使用Bonferroni 校正的成对测试

• 方差分析

当测试两个以上组之间的差异时，我们不得不使用多个测试，这会影响我们的I 型错误率。有几种方法可以控制错误率。我们将看到如何利用 Bonferroni 校正来控制I 型错误率。我们还将在本节中讨论方差分析，它用于测试多个组之间的均值差异。

多次显著性测试

在前几节中，我们探讨了如何比较两组。在本节中，我们将考虑在存在两个以上组时如何进行测试。让我们再次考虑工厂示例，其中我们在工厂地面上有几种机器模型（模型 A、模型 B 和模型 C），这些机器在工厂中用于执行相同的操作。一个有趣的问题可能是：是否有一种机器模型比其他两种模型的平均输出更高？ 为了做出这个判断，我们需要进行三个测试，比较每个模型的均值差异与其他模型，测试均值差异是否不同于零。这些是我们需要测试的零假设：

μ_output,A − μ_output,B = 0

μ_output,B − μ_output,C = 0

μ_output,A − μ_output,C = 0

在进行多重测试时，我们需要应用 p 值校正以控制预期的错误率。回想一下，在 第三章“假设检验”中，我们定义了假设检验的预期错误率为 α。这是我们预期一个 单个 假设检验会导致 I 类 错误的比率。在我们的工厂机器示例中，我们进行了三个假设检验，这意味着 我们犯 I 类错误的概率是三倍。虽然我们的例子具体考虑了均值差异的多重测试，但这适用于任何类型的假设检验。在这些多重测试的情况下，我们通常会定义一个 家族错误率（FWER）并应用 p 值校正以控制 FWER。FWER 是从一组假设检验中犯 I 类 错误的概率。该组内测试的错误率是 个体错误率（IER）。我们将如下定义 IER 和 FWER：

• IER: 个体假设检验的预期 I 类 错误率

• FWER: 一组假设检验的预期 I 类 错误率

我们将在本节中讨论一种 p 值校正方法，以提供对推理的直观理解。

Bonferroni 校正

调整 p 值以控制多重假设检验的一种方法是 Bonferroni 校正。Bonferroni 校正通过均匀降低测试系列中每个单个测试的显著性水平来控制 FWER。鉴于我们在一个系列中有 m 个测试，每个测试的 p 值为 p_i，那么 p 值校正如下所示：

p_i ≤ α_m

以我们之前提到的三种机器模型为例，我们有一个包含三个测试的测试系列，使得 m = 3。如果我们让 FWER 为 0.05，那么，使用 Bonferroni 校正，三个单个测试的显著性水平如下：

0.05 * 3 = 0.0167

因此，在这个例子中，任何单个测试都需要有一个 p 值为 0.0167 才能被认为是显著的。

对 I 类和 II 类错误的影响

如前所述，Bonferroni 校正降低了单个测试的显著性水平，以控制家族级别的 I 类 错误率。我们还应该考虑这种变化如何影响 II 类 错误率。一般来说，降低单个测试的显著性水平会增加该测试犯 II 类 错误的可能性（正如 Bonferroni 校正所做的那样）。虽然我们只在本节中讨论了 Bonferroni 校正，但还有其他 p 值校正方法，它们提供了不同的权衡。请查看 statsmodels 中 multipletests 的文档，以查看 statsmodels 中实现的 p 值校正列表。

让我们通过一个例子来看看如何使用来自UCI 机器学习仓库的Auto MPG数据集的每加仑英里数（MPG）数据，该数据集的链接为：archive.ics.uci.edu/ml/datasets/Auto+MPG [1]。此数据集包含各种属性，包括origin、mpg、cylinders和displacement，涵盖了 1970 年至 1982 年间生产的车辆。我们在这里展示分析的简略形式；完整分析包含在本章代码仓库中的笔记本中。

对于这个例子，我们将使用mpg和origin变量，并在显著性水平 0.01 下测试不同来源的mpg是否存在差异。组均值如下表所示（在此数据集中，origin是一个整数编码标签）。

origin	mpg
1	20.0
2	27.9
3	30.5

图 4.15 – 每个来源组的车辆 MPG 均值

对每个均值进行 t 检验，我们得到以下 p 值：

零假设	未校正的 p 值
μ1 − μ2 = 0	7.946116336281346e-12
μ1 − μ3 = 0	4.608511957238898e-19
μ2 − μ3 = 0	0.0420926104552266

图 4.16 – 每组之间 t 检验的未校正 p 值

对 p 值应用 Bonferroni 校正，我们得到以下 p 值：

零假设	校正后的 p 值（Bonferroni）
μ1 − μ2 = 0	2.38383490e-11
μ1 − μ3 = 0	1.38255359e-18
μ2 − μ3 = 0	0.126277831

图 4.17 – 每组之间 t 检验的校正 p 值

在先前的 p 值基础上，在显著性水平 0.01 下，拒绝组 1 和组 2 均值差异以及组 1 和组 3 均值差异的零假设，但未能拒绝组 2 和组 3 均值差异的零假设。

其他 p 值校正方法

在本节中，我们只讨论了一种 p 值校正方法——Bonferroni 校正，以提供对 p 值校正理由的直观理解。然而，还有其他校正方法可能更适合你的问题。要查看statsmodels中实现的 p 值校正方法列表，请检查statsmodels.stats.multitest.multipletests的文档。

方差分析

在上一节关于多重显著性检验中，我们看到了如何执行多重检验以确定组间均值是否存在差异。当处理组均值时，一个有用的第一步是进行方差分析。ANOVA 是一种用于确定多个组均值之间是否存在差异的统计检验。零假设是没有均值差异，备择假设是均值不全相等。由于 ANOVA 检验均值差异，它通常在执行成对假设检验之前使用。如果 ANOVA 的零假设未能被拒绝，则不需要执行成对检验。然而，如果 ANOVA 的零假设被拒绝，则可以进行成对检验以确定哪些特定的均值存在差异。

ANOVA 与成对检验的比较

虽然成对检验是检验组间差异的一般程序，但 ANOVA 只能用于检验均值差异。

在这个例子中，我们将再次考虑来自 Auto MPG 数据集的车辆 MPG。由于我们已运行成对检验并发现基于来源的车辆平均 mpg 存在显著差异，我们预计 ANOVA 将提供积极的检验结果（拒绝零假设）。执行 ANOVA 计算后，我们得到以下输出。小的 p 值表明我们应该拒绝零假设：

anova = anova_oneway(data.mpg, data.origin, use_var='equal')
print(anova)
# statistic = 98.54179491075868
# pvalue = 1.915486418412936e-35

这里展示的方差分析（ANOVA）是简化的。要查看完整的代码，请参阅本章代码仓库中相关的笔记本。

在本节中，我们介绍了对多于两组数据进行假设检验的方法。第一种方法是带有 p 值校正的成对检验，这是一种通用方法，可用于任何类型的假设检验。我们介绍的另一种方法是 ANOVA，这是一种用于检验组均值差异的特定检验。这不是一种像成对检验那样的通用方法，但可以用作在执行均值差异的成对检验之前的一个步骤。在下一节中，我们将介绍另一种可以用来确定两组数据是否相关的参数检验类型。

皮尔逊相关系数

皮尔逊相关系数，也称为皮尔逊 r（当应用于总体数据时称为皮尔逊 rho (ρ)）或皮尔逊积矩样本相关系数（PPMCC），是一种双变量检验，用于衡量两个变量之间的线性相关性。该系数产生一个介于-1 到 1 之间的值，其中-1 表示强烈的负相关性，1 表示强烈的正相关性。零值的系数表示两个变量之间没有相关性。弱相关性通常被认为是介于 +/- 0.1 和 +/- 0.3 之间，中等相关性是介于 +/- 0.3 和 +/- 0.5 之间，强相关性是介于 +/- 0.5 到 +/- 1.0 之间。

这个测试被认为是参数测试，但不需要假设正态分布或方差齐性。然而，数据必须是独立采样的（既随机选择又无序列相关性），具有有限的方差——例如具有非常重尾部的分布——并且是连续数据类型。该测试不指示输入变量和响应变量；它只是两个变量之间线性关系的度量。该测试使用标准化协方差来推导相关性。回想一下，标准化需要将一个值除以标准差。

这里展示了总体皮尔逊系数ρ的方程：

ρ =  σ xy _ σ x σ y

这里，σ xy 是总体协方差，计算如下：

σ xy =  ∑ i=1 N  (x i − μ x)(y i − μ y)  ______________ N

这里展示了样本皮尔逊系数 r 的方程：

r =  S xy _ S x S y

这里，S xy 是样本协方差，计算如下：

S xy =  ∑ i=1 n  (x i −  _ x )(y i −  _ y )  _____________  n − 1

在 Python 中，我们可以使用scipy的scipy.stats.pearsonr函数执行此测试。在下面的代码片段中，我们使用numpy生成两个正态分布的随机数数据集。我们想要测试两个组之间是否存在相关性，因为它们有一些显著的重叠：

from scipy.stats import pearsonr
import matplotlib.pyplot as plt
import scipy.stats as stats
import seaborn as sns
import pandas as pd
import numpy as np
mu1, sigma1 = 0, 1.1
normally_distributed_1 = np.random.normal(mu1, sigma1, 1000)
mu2, sigma2 = 0, 0.7
normally_distributed_2 = np.random.normal(mu2, sigma2,
    1000)
df_norm = pd.DataFrame({'Distribution':['Distribution 1' for i in range(len(normally_distributed_1))] + ['Distribution 2' for i in range(len(normally_distributed_2))], 'X':np.concatenate([normally_distributed_1, normally_distributed_2])})

在图 4**.18中，我们可以观察到两个相关分布的重叠方差：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_04_018.jpg

图 4.18 – 相关分布

在图4**.18中显示的图中，我们可以看到总体之间的重叠。现在，我们想要使用pearsonr()函数测试相关性，如下所示：

p, r = pearsonr(df_norm.loc[df_norm['Distribution'] == 'Distribution 1', 'X'], df_norm.loc[df_norm['Distribution'] == 'Distribution 2', 'X'])
print("p-value = %.4f"%p)
print("Correlation coefficient = %.4f"%r)

下面的输出表明，在 0.05 的显著性水平下，我们有一个 0.9327 的相关性水平（p 值为 0.0027）：

p-value = 0.0027

相关系数 = 0.9327

为了以不同的方式表达相关性，我们可以说，在分布 2 中，由分布 1 解释的方差水平（r 2，也称为拟合优度或确定系数）为 0.9327 2 = 87%，假设我们知道分布 2 是对分布 1 的反应。否则，我们可以说两个变量之间存在 0.93 的相关性或 87%的方差解释水平。

现在，让我们看看 R 中的Motor Trend Car Road Tests数据集，我们使用statsmodels datasets.get_rdataset函数导入它。这里，我们有前五行，它们包含每加仑英里数（mpg）、气缸数（cyl）、发动机排量（disp）、马力（hp）、后轴齿轮比（drat）、重量（wt）、驾驶四分之一英里所需的最短时间（qsec）、发动机形状（vs=0表示 v 形，vs=1表示直列），变速器（am=0表示自动，am=1表示手动）、齿轮数（gear）和化油器数（carb）（如果不是喷射燃料）：

import statsmodels.api as sm
df_cars = sm.datasets.get_rdataset("mtcars","datasets").data

在图 4.19 中，我们可以看到数据集的前五行，其中包含适合皮尔逊相关分析的数据。

mpg	cyl	disp	hp	drat	wt	qsec	Vs	am	gear	carb
21	6	160	110	3.9	2.62	16.46	0	1	4	4
21	6	160	110	3.9	2.875	17.02	0	1	4	4
22.8	4	108	93	3.85	2.32	18.61	1	1	4	1
21.4	6	258	110	3.08	3.215	19.44	1	0	3	1
18.7	8	360	175	3.15	3.44	17.02	0	0	3	2

图 4.19 – mtcars 数据集的前五行

使用数据集，我们可以使用以下代码绘制相关矩阵，该代码显示了数据集中所有特征的成对相关，以了解它们之间的关系：

sns.set_theme(style="white")
corr = df_cars.corr()
f, ax = plt.subplots(figsize=(15, 10))
cmap = sns.diverging_palette(250, 20, as_cmap=True)
sns.heatmap(corr, cmap=cmap, vmax=.4, center=0,
            square=True, linewidths=.5, annot=True)

假设我们对于对四分之一英里时间（qsec）最有意义的变量感到好奇。在图 4**.20中，通过观察qsec的线条，我们可以看到vs（V 形）与 0.74 的正相关。由于这是一个二元变量，我们可以根据这个数据集假设，直列发动机比 V 形发动机更快。然而，还有其他与显著相关的协变量。例如，与发动机形状几乎同样强烈相关的是马力，即随着马力的增加，四分之一英里运行时间会下降：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/bd-stat-mdl-py/img/B18945_04_020.jpg

图 4.20 – 相关矩阵热图

相关矩阵对于探索多个变量之间的关系非常有用。它也是构建统计和机器学习（ML）模型，如线性回归时进行特征选择的有用工具。

功率分析示例

功率分析是一种统计方法，用于确定假设检验所需的适当样本大小，以便在防止第二类错误（即当零假设应该被拒绝时未能拒绝零假设）时具有足够的功效。功率分析还可以根据样本大小确定可检测的效果量（或差异）之间的差异。换句话说，基于特定的样本大小和分布，功率分析可以为分析师提供研究人员可能能够可靠地识别的特定最小差异。在本节中，我们将通过单样本 t 检验演示功率分析。

单样本 t 检验

假设一家制造商销售一种每月能生产 100,000 个单位的机器，标准差为 2,800 个单位。一家公司购买了这些机器中的一批，并发现它们只能生产 90,000 个单位。该公司想知道需要多少台机器才能以高置信度确定这些机器无法生产 100,000 个单位。以下功率分析表明，对于 t 检验，需要三个机器的样本，以 85% 的概率防止未能识别出实际机器性能与市场宣传机器性能之间的统计显著差异：

from statsmodels.stats.power import TTestPower
import numpy as np
# Difference of distribution mean and the value to be assessed divided by the distribution standard deviation
effect_size = abs(100000-90000) / 2800
powersTT = TTestPower()
result = powersTT.solve_power(effect_size, nobs=3, alpha=0.05, alternative='two-sided')
print('Power based on sample size:{}'.format(round(result,2)))
# Power based on sample size: 0.85

额外的功率分析示例

要查看 Python 中功率分析的更多示例，请参阅本书的 GitHub 仓库。在那里，我们有额外的 t 检验和 F 检验的示例，这些示例专注于分析样本组之间的方差。

摘要

本章涵盖了参数检验的主题。从参数检验的假设开始，我们确定了测试这些假设被违反的方法，并讨论了在所需假设未满足时可以假设稳健性的场景。然后，我们研究了 z 检验的最流行替代方案之一，即 t 检验。我们通过多次应用此测试，包括使用合并、配对和 Welch 的非合并版本进行双样本分析的单样本和双样本版本的测试。接下来，我们探讨了方差分析技术，其中我们研究了如何使用来自多个组的数据来识别它们之间的统计显著差异。这包括当存在大量组时对 p 值的最流行调整之一——Bonferroni 校正，它有助于在执行多个测试时防止 I 类错误 的膨胀。然后，我们研究了使用 Pearson 相关系数对连续数据进行相关性分析，以及如何使用相关性矩阵和相应的热图来可视化相关性。最后，我们简要概述了功率分析，并举例说明了如何使用单样本 t 检验进行此分析。在下一章中，我们将讨论非参数假设检验，包括在本章中与参数检验配对的测试以及当假设不能安全假设时的新测试。

参考文献

[1] Dua, D. 和 Graff, C. (2019). UCI 机器学习仓库 [archive.ics.uci.edu/ml]. Irvine, CA: 加州大学信息与计算机科学学院。