精通 Pandas 第二版（五）
原文：
annas-archive.org/md5/cc15886c82662bbefcd80b2317d3496d译者：飞龙
第十二章：贝叶斯统计与最大似然估计简述
本章将简要介绍一种替代性的统计推断方法，称为贝叶斯统计。这并不是一本完整的入门书籍，而是作为贝叶斯方法的简介。我们还将探讨相关的 Python 库，并学习如何使用pandas和matplotlib来辅助数据分析。以下是将要讨论的各种主题：
-
贝叶斯统计简介
-
贝叶斯统计的数学框架
-
概率分布
-
贝叶斯统计与频率统计
-
PyMC 和蒙特卡洛模拟简介
-
贝叶斯分析示例 – 转折点检测
贝叶斯统计简介
贝叶斯统计学的领域建立在 18 世纪统计学家、哲学家和长老会牧师托马斯·贝叶斯的工作基础上。他著名的贝叶斯定理，构成了贝叶斯统计学的理论基础，并在 1763 年死后发表，作为逆向概率问题的解决方案。更多详细信息请参见en.wikipedia.org/wiki/Thomas_Bayes。
逆向概率问题在 18 世纪初非常流行，通常以以下形式提出。
假设你和朋友玩一个游戏。袋子 1 里有 10 个绿色球和 7 个红色球，袋子 2 里有 4 个绿色球和 7 个红色球。你的朋友抛硬币（没有告诉你结果），随机从一个袋子里抽出一个球，并展示给你。球是红色的。那么，球是从袋子 1 里抽出来的概率是多少？
这些问题被称为逆向概率问题，因为我们试图根据随后的结果（球是红色的）来估计一个已经发生的事件的概率（球是从哪个袋子里抽出来的）：
贝叶斯球示意图
让我们快速说明如何解决之前提到的逆向概率问题。我们想计算在知道球是红色的情况下，球是从袋子 1 里抽出来的概率。这可以表示为P(袋子 1 | 红球)。
让我们从计算选择红球的概率开始。这个概率可以通过按照前述图中的两条红色路径来计算。因此，我们有以下结果：
现在，选择袋子 1 中的红球的概率只通过走上面这条路径来展示，计算方法如下：
从袋子 2 中选择红球的概率计算如下：
注意，这个概率也可以写成如下形式：
我们可以看到，P(袋子 1) = 1/2，树的最后一个分支仅在球首先在袋子 1 中并且是红球时才会被遍历。因此，我们将得到以下结果：
贝叶斯统计的数学框架
贝叶斯方法是一种替代的统计推断方法。我们将首先介绍贝叶斯定理，这是所有贝叶斯推断的基本方程。
在我们开始之前，需要做一些关于概率的定义：
-
A，B：这些是可以以某种概率发生的事件。
-
P(A) 和 P(B)：这是某一特定事件发生的概率。
-
P(A|B)：这是在给定 B 发生的情况下，A 发生的概率。这称为条件概率。
-
P(AB) = P(A 和 B)：这是 A 和 B 同时发生的概率。
我们从以下基本假设开始：
P(AB) = P(B) * P(A|B)
上述方程显示了P(AB)的联合概率与条件概率P(A|B)以及所谓的边际概率P(B)之间的关系。如果我们重写方程，就得到了条件概率的如下表达式：
P(A|B) = P(AB)/P(B)
这在某种程度上是直观的——给定 B 的情况下，A的概率是通过将A和B同时发生的概率除以 B 发生的概率得到的。其思想是，已知 B 发生，因此我们除以其概率。有关该方程更严格的处理，可以参考bit.ly/1bCYXRd，标题为概率：联合、边际和条件概率。
类似地，由于对称性，我们有P(AB) = P(BA) = P(A) * P(B|A)。因此，我们有P(A) * P(B|A) = P(B) * P(A|B)。通过在两边同时除以P(B)并假设P(B) != 0，我们得到如下结果：
上述公式被称为贝叶斯定理，它是所有贝叶斯统计推断的基石。为了将贝叶斯定理与推理统计联系起来，我们将该方程改写为所谓的历时性 解释，如下所示：
这里，H代表假设，D代表已经发生的事件，我们在统计研究中使用它，也称为数据。
表达式 (H) 是我们在观察数据之前对假设的概率。这被称为 先验概率。贝叶斯统计学家常常将先验概率作为一个优势来宣传，因为先前结果的先验知识可以作为当前模型的输入，从而提高准确性。有关更多信息，请参考 www.bayesian-inference.com/advantagesbayesian。
P(D) 是无论假设如何，观察到的数据的概率。这被称为 归一化常数。归一化常数并不总是需要计算，特别是在许多流行的算法中，如 马尔科夫链蒙特卡洛（MCMC），我们将在本章稍后讨论。
P(H|D) 是给定我们观察到的数据后，假设成立的概率。这被称为 后验。
P(D|H) 是考虑我们的假设后，获得数据的概率。这被称为 似然性。
因此，贝叶斯统计就是应用贝叶斯定理来解决推理统计问题，H 代表我们的假设，D 代表数据。
贝叶斯统计模型是以参数为基础的，这些参数的不确定性通过概率分布来表示。这与频率主义方法不同，后者将值视为确定性的。另一种表示方式如下：
这里，Θ 代表我们的未知数据，x 代表我们的观察数据。
在贝叶斯统计中，我们对先验数据做出假设，并利用似然性通过贝叶斯定理更新后验概率。作为说明，下面是一个经典问题，也被称为 urn 问题：
-
两个 urn 中包含有色球
-
urn 1 包含 50 颗红球和 50 颗蓝球
-
urn 2 包含 30 颗红球和 70 颗蓝球
-
从两个 urn 中随机选择一个（50% 的概率），然后从其中一个 urn 中随机抽取一颗球
如果抽到一颗红球，那么它来自 urn 1 的概率是多少？我们想要的是 P(H|D) —— 即 P(球来自 urn 1 | 抽到红球)。
这里，H 表示球是从 urn 1 中抽出的，D 表示抽到的球是红色的：
我们知道 P(H|D) = P(H) * P(D|H)/P(D), P(D|H) = 0.5, P(D) = (50 + 30)/(100 + 100) = 0.4。这也可以表述为：
因此，我们得出结论 P(H|D) = 0.5 * 0.5/0.4 = 0.25/0.4 = 0.625。
贝叶斯理论和赔率
贝叶斯定理有时可以通过使用一种称为赔率的替代概率公式来以更自然和方便的形式表示。赔率通常以比率的形式表示，并且被广泛使用。马匹赢得比赛的 3 比 1 赔率（通常写作 3:1）表示马匹有 75%的概率会获胜。
给定一个概率p，赔率可以计算为 odds = p:(1 - p)，当p=0.75时，结果为 0.75:0.25，即 3:1。
使用赔率，我们可以将贝叶斯定理重写如下：
贝叶斯统计的应用
贝叶斯统计可以应用于我们在经典统计学中遇到的许多问题，如下所示：
-
参数估计
-
预测
-
假设检验
-
线性回归
学习贝叶斯统计学有许多令人信服的理由，比如使用先验信息来更好地优化当前模型。贝叶斯方法基于概率分布而非点估计，因此能够提供更现实的预测。贝叶斯推断基于现有数据对假设进行推理——P(假设|数据)。而频率主义方法则尝试根据假设拟合数据。可以说，贝叶斯方法更具逻辑性和实证性，因为它试图基于事实而非假设来建立信念。有关更多信息，请参见www.bayesian-inference.com/advantagesbayesian。
概率分布
在本节中，我们将简要介绍各种概率分布的特性。许多这些分布用于贝叶斯分析，因此在我们继续之前，需要简要概述它们。我们还将演示如何使用matplotlib生成和显示这些分布。为了避免在每个代码片段中重复import语句，我们将在以下命令中展示一组标准的 Python 代码导入，在任何后续代码片段中都需要先运行这些导入。每个会话只需运行一次这些导入。导入语句如下：
import pandas as pd
import numpy as np
from matplotlib import pyplot as plt
from matplotlib import colors
import matplotlib.pyplot as plt
import matplotlib
%matplotlib inline
拟合分布
在贝叶斯分析中，我们必须进行的一个步骤是将我们的数据拟合到一个概率分布中。选择正确的分布有时需要一些技巧，通常需要统计知识和经验，但我们可以遵循一些指导原则来帮助我们。以下是这些指导原则：
-
确定数据是离散的还是连续的
-
检查数据的偏斜/对称性，如果有偏斜，确定其方向
-
确定是否存在上下限
-
确定在分布中观察到极端值的可能性
统计试验是一个可重复的实验，具有一组已定义的结果，这些结果组成了样本空间。伯努利试验是一个“是/否”实验，其中随机X变量在“是”的情况下赋值为 1，在“否”的情况下赋值为 0。抛硬币并查看是否正面朝上的实验就是一个伯努利试验的例子。
概率分布有两种类型：离散分布和连续分布。在接下来的部分中，我们将讨论这两类分布的区别，并了解主要的分布。
离散概率分布
在这种情况下，变量只能取特定的离散值，如整数。一个离散随机变量的例子是，当我们抛掷硬币五次时获得的正面次数：可能的值为{0,1,2,3,4,5}——例如，不能获得 3.82 次正面。随机变量可以取的值的范围由被称为概率质量函数（PMF）的函数来指定。
离散均匀分布
离散均匀分布是一种分布，用于模拟具有有限可能结果的事件，其中每个结果的出现概率相等。对于 n个结果，每个结果的发生概率为1/n。
这方面的一个例子是投掷一个公平的骰子。六个结果中任何一个的概率都是1**/6。概率质量函数（PMF）由1/n给出，期望值和方差分别由 (max + min)/2 和 (n²-1)/12 给出：
from matplotlib import pyplot as plt
import matplotlib.pyplot as plt
X=range(0,11)
Y=[1/6.0 if x in range(1,7) else 0.0 for x in X]
plt.plot(X,Y,'go-', linewidth=0, drawstyle='steps-pre',
label="p(x)=1/6")
plt.legend(loc="upper left")
plt.vlines(range(1,7),0,max(Y), linestyle='-')
plt.xlabel('x')
plt.ylabel('p(x)')
plt.ylim(0,0.5)
plt.xlim(0,10)
plt.title('Discrete uniform probability distribution with
p=1/6')
plt.show()
以下是输出：
离散均匀分布
伯努利分布
伯努利分布用于衡量试验中的成功概率，例如，抛硬币正面或反面朝上的概率。它可以通过一个随机X变量来表示，如果硬币是正面朝上，则取值为 1；如果是反面朝上，则取值为 0。正面朝上的概率用p表示，反面朝上的概率用q=1-p表示。
这可以通过以下概率质量函数（PMF）表示：
期望值和方差由以下公式给出：
如需更多信息，请访问 en.wikipedia.org/wiki/Bernoulli_distribution。
我们现在将使用matplotlib和scipy.stats绘制伯努利分布，如下所示：
In [20]:import matplotlib
from scipy.stats import bernoulli
a = np.arange(2)
colors = matplotlib.rcParams['axes.color_cycle']
plt.figure(figsize=(12,8))
for i, p in enumerate([0.0, 0.2, 0.5, 0.75, 1.0]):
ax = plt.subplot(1, 5, i+1)
plt.bar(a, bernoulli.pmf(a, p), label=p, color=colors[i], alpha=0.5)
ax.xaxis.set_ticks(a)
plt.legend(loc=0)
if i == 0:
plt.ylabel("PDF at $k$")
plt.suptitle("Bernoulli probability for various values of $p$")
以下是输出：
伯努利分布输出
二项分布
二项分布用于表示在n次独立伯努利试验中成功的次数，表达式如下：
Y = X[1] + X[2] + ..**. + X[n]
使用抛硬币的例子，这个分布描述了在n次试验中获得X个正面的概率。对于 100 次抛掷，二项分布模型描述了获得 0 个正面（极不可能）到 50 个正面（最有可能）再到 100 个正面（也极不可能）的几率。这导致当概率完全均等时，二项分布是对称的，而当概率不均等时，则呈偏斜分布。PMF 由以下表达式给出：
期望和方差分别由以下表达式给出：
这通过以下代码块显示：
from scipy.stats import binom
clrs = ['blue','green','red','cyan','magenta'] plt.figure(figsize=(12,6))
k = np.arange(0, 22)
for p, color in zip([0.001, 0.1, 0.3, 0.6, 0.999], clrs):
rv = binom(20, p)
plt.plot(k, rv.pmf(k), lw=2, color=color, label="$p$=" + str(round(p,1)))
plt.legend()
plt.title("Binomial distribution PMF")
plt.tight_layout()
plt.ylabel("PDF at $k$")
plt.xlabel("$k$")
以下是输出结果：
二项分布
泊松分布
泊松分布模型描述了在给定时间间隔内事件发生的概率，假设这些事件以已知的平均速率发生，并且每次事件发生与上次事件发生后的时间间隔无关。
一个可以通过泊松分布建模的具体过程的例子是，如果一个人每天平均收到 23 封电子邮件。假设这些邮件的到达时间是相互独立的，那么一个人每天收到的邮件总数就可以用泊松分布来建模。
另一个例子可以是每小时在特定车站停靠的火车数量。泊松分布的概率质量函数（PMF）由以下表达式给出：
其中，λ是速率参数，表示单位时间内期望的事件/到达次数，而k是表示事件/到达次数的随机变量。
期望和方差分别由以下表达式给出：
更多信息请参考en.wikipedia.org/wiki/Poisson_process。
使用matplotlib绘制不同值的 PMF，如下所示：
In [11]: %matplotlib inline
import numpy as np
import matplotlib
import matplotlib.pyplot as plt
from scipy.stats import poisson
colors = matplotlib.rcParams['axes.color_cycle']
k=np.arange(15)
plt.figure(figsize=(12,8))
for i, lambda_ in enumerate([1,2,4,6]):
plt.plot(k, poisson.pmf(k, lambda_), '-o',
label="$\lambda$=" + str(lambda_), color=colors[i])
plt.legend()
plt.title("Possion distribution PMF for various $\lambda$")
plt.ylabel("PMF at $k$")
plt.xlabel("$k$")
plt.show()
以下是输出结果：
泊松分布
几何分布
几何分布用于独立的伯努利试验，并衡量获得一次成功所需的试验次数（X）。它也可以表示第一次成功之前的失败次数（Y = X-1）。
PMF 表示为以下形式：
上述表达式是有道理的，因为f(k) = P(X =k)，如果成功需要进行k次试验（p），那么这意味着我们必须有k-1次失败，失败的概率是*(1-p)*。
期望和方差分别给出如下：
以下命令清晰地解释了前面的公式：
In [12]: from scipy.stats import geom
p_vals=[0.01,0.2,0.5,0.8,0.9]
x = np.arange(geom.ppf(0.01,p),geom.ppf(0.99,p))
colors = matplotlib.rcParams['axes.color_cycle']
for p,color in zip(p_vals,colors):
x = np.arange(geom.ppf(0.01,p),geom.ppf(0.99,p))
plt.plot(x,geom.pmf(x,p),'-o',ms=8,label='$p$=' + str(p))
plt.legend(loc='best')
plt.ylim(-0.5,1.5)
plt.xlim(0,7.5)
plt.ylabel("Pmf at $k$")
plt.xlabel("$k$")
plt.title("Geometric distribution PMF")
以下是输出结果：
几何分布
负二项分布
负二项分布用于独立的伯努利试验，衡量在发生指定次数的成功（r）之前所需的试验次数（X=k）。一个例子是，投掷硬币直到得到五次正面所需的次数。其概率质量函数（PMF）如下所示：
期望值和方差分别由以下公式给出：
我们可以看到，负二项分布是几何分布的推广，几何分布是负二项分布的一个特例，其中r=1。
代码和图表如下所示：
In [189]: from scipy.stats import nbinom
from matplotlib import colors
clrs = matplotlib.rcParams['axes.color_cycle']
x = np.arange(0,11)
n_vals = [0.1,1,3,6]
p=0.5
for n, clr in zip(n_vals, clrs):
rv = nbinom(n,p)
plt.plot(x,rv.pmf(x), label="$n$=" + str(n), color=clr)
plt.legend()
plt.title("Negative Binomial Distribution PMF")
plt.ylabel("PMF at $x$")
plt.xlabel("$x$")
以下是输出结果：
负二项分布
连续概率分布
在连续概率分布中，变量可以取任何实数值。它不像离散概率分布那样被限制为有限的取值集合；例如，健康新生儿的平均体重大约在 6 到 9 磅之间（例如，其体重可以是 7.3 磅）。连续概率分布的特点是概率密度函数（PDF）。
随机变量可能取的所有概率之和为 1。因此，概率密度函数图形下的面积加起来为 1。
连续均匀分布
均匀分布模型表示一个随机变量X，其值可以在区间*[a, b]*内均匀分布。
概率密度函数由以下公式给出：，当a ≤ x ≤ b时为非零值，否则为零。
期望值和方差分别由以下公式给出：
以下代码生成并绘制了不同样本大小的连续均匀概率分布：
In [11]: np.random.seed(100) # seed the random number generator
# so plots are reproducible
subplots = [111,211,311]
ctr = 0
fig, ax = plt.subplots(len(subplots), figsize=(10,12))
nsteps=10
for i in range(0,3):
cud = np.random.uniform(0,1,nsteps) # generate distrib
count, bins, ignored = ax[ctr].hist(cud,15,normed=True)
ax[ctr].plot(bins,np.ones_like(bins),linewidth=2, color='r')
ax[ctr].set_title('sample size=%s' % nsteps)
ctr += 1
nsteps *= 100
fig.subplots_adjust(hspace=0.4)
plt.suptitle("Continuous Uniform probability distributions for various sample sizes" , fontsize=14)
以下是输出结果：
连续均匀分布
指数分布
指数分布描述了在泊松过程中两个事件之间的等待时间。泊松过程是一个符合泊松分布的过程，事件以已知的平均速率不可预测地发生。指数分布可以视为几何分布的连续极限，并且也是马尔科夫过程（无记忆过程）。
一个无记忆随机变量具有以下特性：它的未来状态仅依赖于当前时刻的相关信息，而不依赖于更久远过去的信息。对马尔科夫/无记忆随机变量建模的一个例子是基于随机游走理论对短期股价行为的建模。这引出了金融学中的有效市场假说。更多信息，请参阅en.wikipedia.org/wiki/Random_walk_hypothesis。
指数分布的概率密度函数（PDF）为 f(x) = λe^(-λx)。期望值和方差分别由以下公式给出：
E(X) = 1/λ**Var(X) = 1/λ²
参考资料，请参阅en.wikipedia.org/wiki/Exponential_distribution。
分布图和代码如下所示：
In [15]: import scipy.stats
clrs = colors.cnames
x = np.linspace(0,4, 100)
expo = scipy.stats.expon
lambda_ = [0.5, 1, 2, 5]
plt.figure(figsize=(12,4))
for l,c in zip(lambda_,clrs):
plt.plot(x, expo.pdf(x, scale=1./l), lw=2,
color=c, label = "$\lambda = %.1f$"%l)
plt.legend()
plt.ylabel("PDF at $x$")
plt.xlabel("$x$")
plt.title("Pdf of an Exponential random variable for various $\lambda$");
以下是输出结果：
指数分布
正态分布
在统计学中，最重要的分布可能就是正态/高斯分布。它模拟了围绕一个中心值的概率分布，且没有左右偏差。有许多现象符合正态分布，例如以下这些：
-
婴儿的出生体重
-
测量误差
-
血压
-
测试分数
正态分布的重要性通过中心极限定理得到了强调，该定理指出，从同一分布中独立抽取的许多随机变量的均值大致呈正态分布，无论原始分布的形态如何。其期望值和方差分别如下：
E(X) = μ**Var(X) = σ²
正态分布的概率密度函数（PDF）由以下公式给出：
以下代码和图解释了这个公式：
In [54]: import matplotlib
from scipy.stats import norm
X = 2.5
dx = 0.1
R = np.arange(-X,X+dx,dx)
L = list()
sdL = (0.5,1,2,3)
for sd in sdL:
f = norm.pdf
L.append([f(x,loc=0,scale=sd) for x in R])
colors = matplotlib.rcParams['axes.color_cycle']
for sd,c,P in zip(sdL,colors,L):
plt.plot(R,P,zorder=1,lw=1.5,color=c,
label="$\sigma$=" + str(sd))
plt.legend()
ax = plt.axes()
ax.set_xlim(-2.1,2.1)
ax.set_ylim(0,1.0)
plt.title("Normal distribution Pdf")
plt.ylabel("PDF at $\mu$=0, $\sigma$")
以下是输出结果：
正态分布
用于绘制分布的 Python 代码参考可以在 bit.ly/1E17nYx 找到。
正态分布也可以看作是二项分布的连续极限，以及其他分布，如 。我们可以在以下代码和图中看到这一点：
In [18]:from scipy.stats import binom
from matplotlib import colors
cols = colors.cnames
n_values = [1, 5,10, 30, 100]
subplots = [111+100*x for x in range(0,len(n_values))]
ctr = 0
fig, ax = plt.subplots(len(subplots), figsize=(6,12))
k = np.arange(0, 200)
p=0.5
for n, color in zip(n_values, cols):
k=np.arange(0,n+1)
rv = binom(n, p)
ax[ctr].plot(k, rv.pmf(k), lw=2, color=color)
ax[ctr].set_title("$n$=" + str(n))
ctr += 1
fig.subplots_adjust(hspace=0.5)
plt.suptitle("Binomial distribution PMF (p=0.5) for various values of n", fontsize=14)
以下是输出结果：
随着 n 增加，二项分布趋近于正态分布。事实上，这在前面的 n>=30 的图中可以清楚地看到。
贝叶斯统计与频率派统计
在当今的统计学中，对于如何解释数据并进行统计推断，有两种学派。迄今为止，经典且占主导地位的方法是所谓的频率学派方法（参见第七章，统计学概览 – 经典方法）。我们在这一章中探讨的是贝叶斯方法。
什么是概率？
贝叶斯学派与频率学派世界观之间争论的核心问题是如何定义概率。
在频率学派的世界观中，概率是通过重复事件的频率推导出来的概念——例如，我们将公平硬币投掷时出现正面的概率定义为一半。这是因为当我们反复投掷公平硬币时，正面朝上的次数与总投掷次数的比值会随着投掷次数的增加接近 0.5。
贝叶斯的世界观不同，概率的概念是与一个人对事件发生的信念程度相关。因此，对于贝叶斯统计学家来说，相信公平骰子出现五点的概率是1/6，这与我们对该事件发生的概率信念有关。
模型如何定义
从模型定义的角度来看，频率学派通过反复实验来分析数据和计算指标如何变化，同时保持模型参数不变。贝叶斯学派则利用固定的实验数据，但在模型参数上变化他们的信念程度。具体解释如下：
-
频率学派：如果模型是固定的，数据是变化的。
-
贝叶斯学派：如果数据是固定的，模型是变化的。
频率学派的方法是使用最大似然法来估计模型参数。它涉及从一组独立同分布的观测数据中生成数据，并将观察到的数据拟合到模型中。最能拟合数据的模型参数值就是最大似然估计器（MLE），它有时可以是观测数据的一个函数。
贝叶斯主义从概率框架的角度以不同的方式处理问题。它使用概率分布来描述值的不确定性。贝叶斯实践者通过观察数据来估计概率。为了计算这些概率，他们使用一个单一的估计器，即贝叶斯公式。这产生了一个分布，而不仅仅是像频率学派方法那样的点估计。
置信区间（频率学派）与可信区间（贝叶斯学派）
让我们比较一下 95% 置信区间（频率学派使用的术语）和 95% 可信区间（贝叶斯实践者使用的术语）的含义。
在频率主义框架中，95%的置信区间意味着，如果你将实验重复无限次，每次生成区间，95%的这些区间将包含我们试图估计的参数，通常称为θ。在这种情况下，区间是随机变量，而不是参数估计值θ，后者在频率主义视角中是固定的。
在贝叶斯可信区间的情况下，我们有一个与频率主义置信区间的传统解释更一致的解释。因此，我们可以得出结论，Pr(a(Y) < θ < b(Y)|θ) = 0.95。在这种情况下，我们可以合理地得出结论，θ有 95%的概率位于该区间内。
欲了解更多信息，请参考*《频率主义与贝叶斯主义：有什么大不了的？》*（Jake VanderPlas，SciPy，2014）在www.youtube.com/watch?v=KhAUfqhLakw。
进行贝叶斯统计分析
进行贝叶斯统计分析包括以下步骤：
-
指定概率模型：在这一步中，我们使用概率分布完全描述模型。基于我们所取样本的分布，我们尝试拟合一个模型并为未知参数分配概率。
-
计算后验分布：后验分布是我们基于观测数据计算出的分布。在这种情况下，我们将直接应用贝叶斯公式。它将作为我们在前一步中指定的概率模型的函数进行指定。
-
检查我们的模型：这是一个必要的步骤，我们在做出推断之前，检查我们的模型及其输出。贝叶斯推断方法使用概率分布为可能的结果分配概率。
似然函数的蒙特卡洛估计与 PyMC
贝叶斯统计不仅仅是另一种方法。它是一个完全不同的统计学实践范式。它使用概率模型来进行推断，基于收集到的数据。这可以通过一个基本的表达式来表示，即P(H|D)。
在这里，H是我们的假设，即我们试图证明的内容，D是我们的数据或观察结果。
为了提醒我们之前的讨论，贝叶斯定理的历时形式如下：
在这里，*P(H)*是无条件的先验概率，表示我们在进行实验之前所知道的内容。*P(D|H)*是我们的似然函数，即在假设为真时，获得我们观察到的数据的概率。
P(D)是数据的概率，也称为归一化常数。它可以通过对H上的分子进行积分得到。
似然函数是我们贝叶斯计算中最重要的部分，它封装了关于数据中未知信息的所有信息。它与反向概率质量函数有一些相似之处。
反对采用贝叶斯方法的一个论点是先验的计算可能具有主观性。支持这种方法的有许多论点，其中之一是可以加入外部先验信息，如前所述。
似然值表示一个未知的积分，在简单的情况下可以通过解析积分获得。
蒙特卡洛 (MC) 积分在涉及更复杂的高维积分的使用案例中是必需的，并且可以用于计算似然函数。
MC 积分可以通过多种抽样方法计算，例如均匀抽样、分层抽样和重要性抽样。在蒙特卡洛积分中，我们可以如下近似积分：
以下是有限和：
这里，x 是来自 g 的样本向量。可以通过大数法则并确保仿真误差较小来证明这个估计是有效的。
在 Python 中进行贝叶斯分析时，我们需要一个模块来帮助我们使用蒙特卡洛方法计算似然函数。PyMC 库满足了这个需求。它提供了一种常用的蒙特卡洛方法，称为 MCMC。我们不会进一步探讨 MCMC 的技术细节，但有兴趣的读者可以通过以下来源了解有关 PyMC 中 MCMC 实现的更多信息：
-
贝叶斯估计中的蒙特卡洛积分：
bit.ly/1bMALeu -
马尔科夫链蒙特卡洛最大似然：
bit.ly/1KBP8hH -
使用 Python 进行贝叶斯统计分析–第一部分，SciPy 2014，Chris Fonnesbeck：
www.youtube.com/watch?v=vOBB_ycQ0RA
MCMC 并非万能的灵丹妙药；这种方法存在一些缺点，其中之一就是算法收敛较慢。
贝叶斯分析示例 – 转折点检测
在这里，我们将尝试使用贝叶斯推理并对一个有趣的数据集进行建模。该数据集由作者的Facebook (FB) 发帖历史构成。我们已经清洗了 FB 历史数据并将日期保存在 fb_post_dates.txt 文件中。文件中的数据如下所示：
head -2 ../fb_post_dates.txt
Tuesday, September 30, 2014 | 2:43am EDT
Tuesday, September 30, 2014 | 2:22am EDT
因此，我们看到一个日期时间序列，表示作者在 Facebook 上发帖的日期和时间。首先，我们将文件读取到 DataFrame 中，并将时间戳拆分成 Date 和 Time 两列：
In [91]: filePath="./data/fb_post_dates.txt"
fbdata_df=pd.read_csv(filePath, sep='|', parse_dates=[0], header=None,names=['Date','Time'])
接下来，我们按如下方式检查数据：
In [92]: fbdata_df.head() #inspect the data
Out[92]: Date Time
0 2014-09-30 2:43am EDT
1 2014-09-30 2:22am EDT
2 2014-09-30 2:06am EDT
3 2014-09-30 1:07am EDT
4 2014-09-28 9:16pm EDT
现在，我们按 Date 索引数据，创建一个 DatetimeIndex，以便我们可以对其进行重采样并按月计数，如下所示：
In [115]: fbdata_df_ind=fbdata_df.set_index('Date')
fbdata_df_ind.head(5)
Out[115]:
Date Time
2014-09-30 2:43am EDT
2014-09-30 2:22am EDT
2014-09-30 2:06am EDT
2014-09-30 1:07am EDT
2014-09-28 9:16pm EDT
然后，我们按如下方式显示关于索引的信息：
In [116]: fbdata_df_ind.index
Out[116]: <class 'pandas.tseries.index.DatetimeIndex'>
[2014-09-30, ..., 2007-04-16]
Length: 7713, Freq: None, Timezone: None
接下来，我们使用重采样来按月获取帖子数量：
In [99]: fb_mth_count_=fbdata_df_ind.resample('M', how='count')
fb_mth_count_.rename(columns={'Time':'Count'},
inplace=True) # Rename
fb_mth_count_.head()
Out[99]: Count
Date
2007-04-30 1
2007-05-31 0
2007-06-30 5
2007-07-31 50
2007-08-31 24
Date 格式显示的是每月的最后一天。现在，我们创建一个从 2007 年到 2015 年的 FB 帖子数量散点图，并使点的大小与 matplotlib 中的值成正比：
In [108]: %matplotlib inline
import datetime as dt
#Obtain the count data from the DataFrame as a dictionary
year_month_count=fb_bymth_count.to_dict()['Count']
size=len(year_month_count.keys())
#get dates as list of strings
xdates=[dt.datetime.strptime(str(yyyymm),'%Y%m')
for yyyymm in year_month_count.keys()]
counts=year_month_count.values()
plt.scatter(xdates,counts,s=counts)
plt.xlabel('Year')
plt.ylabel('Number of Facebook posts')
plt.show()
以下是输出结果：
我们希望调查的问题是，是否在某个时刻行为发生了变化。具体来说，我们想要识别是否存在一个特定时期，FB 帖子的平均数量发生了变化。这通常被称为时间序列中的切换点或变化点。
我们可以利用泊松分布来建模这个。你可能记得，泊松分布可以用来建模时间序列的计数数据。（更多信息请参考 bit.ly/1JniIqy。）
如果我们用 C[i] 表示每月的 FB 帖子数，我们可以将我们的模型表示如下：
r[i] 参数是泊松分布的速率参数，但我们不知道它的值。如果我们查看 FB 时间序列的散点图，可以看到 2010 年中到晚些时候，帖子数量出现了跳跃，可能与 2010 年南非世界杯的开始时间重合，而作者也参加了这场世界杯。
s 参数是切换点，即速率参数发生变化的时刻，而 e 和 l 分别是切换点前后 r[i] 参数的值。可以表示如下：
请注意，这里指定的变量—C, s, e, r 和 l—都是贝叶斯随机变量。对于表示某人对其值的信念的贝叶斯随机变量，我们需要使用概率分布来建模它们。我们希望推断 e 和 l 的值，这些值是未知的。在 PyMC 中，我们可以使用随机和确定性类来表示随机变量。我们注意到，指数分布表示的是泊松事件之间的时间。因此，对于 e 和 l，我们选择使用指数分布来建模它们，因为它们可以是任何正数：
对于 s，我们选择使用均匀分布来建模它，这反映了我们认为切换点可能在整个时间段内的任何一天发生的概率是相等的。这意味着我们有如下假设：
在这里，*t[0], t[f]*分别对应年份的下限和上限，i。现在让我们使用PyMC来表示我们之前开发的模型。我们将使用PyMC来看看是否能在 FB 发布数据中检测到转折点。除了散点图，我们还可以在条形图中显示数据。为此，我们首先需要获取按月份排序的 FB 发布计数列表：
In [69]: fb_activity_data = [year_month_count[k] for k in
sorted(year_month_count.keys())]
fb_activity_data[:5]
Out[70]: [1, 0, 5, 50, 24]
In [71]: fb_post_count=len(fb_activity_data)
我们使用matplotlib来渲染条形图：
In [72]: from IPython.core.pylabtools import figsize
import matplotlib.pyplot as plt
figsize(8, 5)
plt.bar(np.arange(fb_post_count),
fb_activity_data, color="#49a178")
plt.xlabel("Time (months)")
plt.ylabel("Number of FB posts")
plt.title("Monthly Facebook posts over time")
plt.xlim(0,fb_post_count);
以下是输出结果：
从前面的条形图来看，我们能否得出结论认为 FB 频率发布行为在一段时间内发生了变化？我们可以在已开发的模型上使用PyMC来帮助我们找出变化，方法如下：
In [88]: # Define data and stochastics
import pymc as pm
switchpoint = pm.DiscreteUniform('switchpoint',
lower=0,
upper=len(fb_activity_data)-1,
doc='Switchpoint[month]')
avg = np.mean(fb_activity_data)
early_mean = pm.Exponential('early_mean', beta=1./avg)
late_mean = pm.Exponential('late_mean', beta=1./avg)
late_mean
Out[88]:<pymc.distributions.Exponential 'late_mean' at 0x10ee56d50>
在这里，我们为速率参数r定义了一个方法，并使用之前讨论的泊松分布来建模计数数据：
In [89]: @pm.deterministic(plot=False)
def rate(s=switchpoint, e=early_mean, l=late_mean):
''' Concatenate Poisson means '''
out = np.zeros(len(fb_activity_data))
out[:s] = e
out[s:] = l
return out
fb_activity = pm.Poisson('fb_activity', mu=rate,
value=fb_activity_data, observed=True)
fb_activity
Out[89]: <pymc.distributions.Poisson 'fb_activity' at 0x10ed1ee50>
在前面的代码片段中，@pm.deterministic是一个装饰器，表示速率函数是确定性的，即其值完全由其他变量（在本例中为e、s和l）决定。该装饰器是必要的，它告诉PyMC将速率函数转换为确定性对象。如果我们不指定装饰器，将会发生错误。有关 Python 装饰器的更多信息，请参考bit.ly/1zj8U0o。
更多信息，请参考以下网页：
现在，我们将使用 FB 计数数据（fb_activity）和e, s, l（分别为early_mean、late_mean和rate）参数来创建模型。
接下来，使用PyMC，我们创建一个MCMC对象，使我们能够使用 MCMC 方法来拟合数据。然后，我们调用该MCMC对象上的 sample 函数进行拟合：
In [94]: fb_activity_model=pm.Model([fb_activity,early_mean,
late_mean,rate])
In [95]: from pymc import MCMC
fbM=MCMC(fb_activity_model)
In [96]: fbM.sample(iter=40000,burn=1000, thin=20)
[-----------------100%-----------------] 40000 of 40000
complete in 11.0 sec
使用 MCMC 拟合模型涉及利用马尔可夫链蒙特卡洛方法生成后验的概率分布，P(s,e,l | D)。它使用蒙特卡洛过程反复模拟数据采样，并在算法似乎收敛到稳定状态时停止，收敛依据多个标准。这是一个马尔可夫过程，因为连续的样本仅依赖于前一个样本。有关马尔可夫链收敛的更多信息，请参考bit.ly/1IETkhC。
生成的样本被称为轨迹。我们可以通过查看其轨迹的直方图来观察参数的边际后验分布的形态：
In [97]: from pylab import hist,show
%matplotlib inline
hist(fbM.trace('late_mean')[:])
Out[97]: (array([ 15., 61., 214., 421., 517., 426., 202.,
70., 21., 3.]),
array([ 102.29451192, 103.25158404, 104.20865616,
105.16572829, 106.12280041, 107.07987253,
108.03694465, 108.99401677, 109.95108889,
110.90816101, 111.86523313]),
<a list of 10 Patch objects>)
以下是输出结果：
接下来，我们计算早期的均值：
In [98]:plt.hist(fbM.trace('early_mean')[:]) Out[98]: (array([ 20., 105., 330., 489., 470., 314., 147.,
60., 3., 12.]),
array([ 49.19781192, 50.07760882, 50.95740571,
51.83720261, 52.71699951, 53.59679641,
54.47659331, 55.35639021, 56.2361871 ,
57.115984 , 57.9957809 ]),
<a list of 10 Patch objects>)
以下是输出结果：
在这里，我们看到切换点在月份数上的表现：
In [99]: fbM.trace('switchpoint')[:] Out[99]: array([38, 38, 38, ..., 35, 35, 35])
In [150]: plt.hist(fbM.trace('switchpoint')[:]) Out[150]: (array([ 1899., 0., 0., 0., 0., 0.,
0., 0., 0., 51.]),
array([ 35\. , 35.3, 35.6, 35.9, 36.2, 36.5, 36.8,
37.1, 37.4, 37.7, 38\. ]),
<a list of 10 Patch objects>)
以下是输出结果：
切换点的月份数直方图
我们可以看到，切换点大致位于 35 到 38 个月之间。这里，我们使用matplotlib来展示e、s和l的边际后验分布，并将它们绘制在一张图中：
In [141]: early_mean_samples=fbM.trace('early_mean')[:]
late_mean_samples=fbM.trace('late_mean')[:]
switchpoint_samples=fbM.trace('switchpoint')[:]
In [142]: from IPython.core.pylabtools import figsize
figsize(12.5, 10)
# histogram of the samples:
fig = plt.figure()
fig.subplots_adjust(bottom=-0.05)
n_mths=len(fb_activity_data)
ax = plt.subplot(311)
ax.set_autoscaley_on(False)
plt.hist(early_mean_samples, histtype='stepfilled',
bins=30, alpha=0.85, label="posterior of $e$",
color="turquoise", normed=True)
plt.legend(loc="upper left")
plt.title(r"""Posterior distributions of the variables
$e, l, s$""",fontsize=16)
plt.xlim([40, 120])
plt.ylim([0, 0.6])
plt.xlabel("$e$ value",fontsize=14)
ax = plt.subplot(312)
ax.set_autoscaley_on(False)
plt.hist(late_mean_samples, histtype='stepfilled',
bins=30, alpha=0.85, label="posterior of $l$",
color="purple", normed=True)
plt.legend(loc="upper left")
plt.xlim([40, 120])
plt.ylim([0, 0.6])
plt.xlabel("$l$ value",fontsize=14)
plt.subplot(313)
w = 1.0 / switchpoint_samples.shape[0] *
np.ones_like(switchpoint_samples)
plt.hist(switchpoint_samples, bins=range(0,n_mths), alpha=1,
label=r"posterior of $s$", color="green",
weights=w, rwidth=2.)
plt.xlim([20, n_mths - 20])
plt.xlabel(r"$s$ (in days)",fontsize=14)
plt.ylabel("probability")
plt.legend(loc="upper left")
plt.show()
以下是输出结果：
边际后验分布
PyMC还具有绘图功能，因为它使用了matplotlib。在以下图表中，我们展示了时间序列图、自相关图（acorr）和为早期均值、晚期均值以及切换点所绘制的样本的直方图。直方图有助于可视化后验分布。自相关图显示了前一时期的值是否与当前时期的值有较强的关系：
In [100]: from pymc.Matplot import plot
plot(fbM)
Plotting late_mean
Plotting switchpoint
Plotting early_mean
以下是晚期均值图：
pymc_comprehensive_late_mean 的图表
在这里，我们展示了切换点图：
PyMC 综合切换点
在这里，我们展示了早期均值图：
PyMC 综合早期均值
从PyMC的输出结果中，我们可以得出结论，切换点大约在时间序列开始后的 35 到 38 个月之间。这个时间段大致对应于 2010 年 3 月到 7 月之间。作者可以证实，这是他使用 Facebook 的一个标志性年份，因为那年是南非举办 FIFA 世界杯决赛的年份，而他也亲自参加了。
最大似然估计
最大似然估计（MLE）是一种用于从现有样本数据中估计总体分布参数的方法。MLE 方法也可以被视为贝叶斯的替代方法。
概率分布给出了在给定分布参数（如均值、标准差和自由度）下观察到数据点的概率。
给定分布参数下数据点的概率表示为 Prob(X|µ,α) -------1。
MLE 处理的是逆问题。它用于在给定数据点的情况下，找到最可能的分布参数值。为此，定义了一个称为似然度的统计量。似然度被定义为在给定数据点的情况下观察到分布参数的概率。
给定数据点的分布参数的概率表示为 L(µ,α|X)----2。
方程 1 和 2 中的量是相同的概率，只是表述方式不同。因此，我们可以写出以下公式：
为了更好地理解这个概念，请看以下图表：
正态分布 MLE 估计示例，包含两个数据点
该图展示了三个具有不同均值和标准差的正态分布。两条竖直线分别表示值 V[1]=-2.5 和 V[2]=2.5。
假设我们表示数据点 V[1]=-2.5 属于红色分布（均值=0，标准差=1.5）的概率为P(red|V[1]=-2.5)。类似地，数据点 V[1]=-2.5 属于蓝色分布（均值=0，标准差=4）的概率为P(blue|V[1]=-2.5)。
现在，查看这里展示的图表，我们可以得出以下结论：
如果我们必须仅根据现有数据点做出决策，那么我们会决定V[1]=2.5 属于蓝色分布，因为V[1]B[1] > V[1]R[1] > V[1]G[1]，我们选择具有该数据点属于该分布最大概率的分布。
但是如果我们有一个或更多的数据点呢？对于这种情况，让我们向数据集中添加另一个数据点，称为V[2]。V[2]的单个概率如下：
由于现在有两个数据点，我们不能仅凭单个概率做出决策，而必须计算联合概率。如果我们假设一个数据点的发生事件与另一个数据点的发生事件是独立的，那么联合概率将等于单个概率。
假设给定两个数据点时，这两个点属于红色分布的联合概率表示为P(red|V[1]=-2.5, V[2]=2.5)。那么以下等式成立：
我们应该选择具有最大联合概率的分布。在这种情况下，蓝色的分布具有最高的联合概率，我们可以从图表中得出这一结论。
随着数据集中的点数增加，单凭观察前面的图表已经无法再有效地推断最大联合概率。我们需要求助于代数和微积分方法，找出能够最大化数据集属于某个分布的联合概率的分布参数。
如果我们假设单个数据点之间是独立的，那么在给定所有数据点的情况下观察分布参数的似然性或概率可以通过以下公式表示：
在最大似然估计（MLE）计算中，我们尝试找到能够最大化 L(µ,α|X[1], X[2], X[3], ...., X[n]) 的 µ 和 α 的值。在这个过程中，对两边取对数非常有帮助。因为对数是单调递增的函数，它不会改变目标函数，但能使计算变得更加简单。似然函数的对数通常被称为对数似然函数，计算方式如下：
为了找到对数似然函数的最大值，log(L(µ,α|X[1], X[2], X[3], ...., X[n]))，我们可以执行以下操作：
-
对 log(L(µ,α|X[1], X[2], X[3], ...., X[n])) 函数关于 µ 和 α 求一阶导数，并令其等于零
-
对 log(L(µ,α|X[1], X[2], X[3], ...., X[n])) 函数关于 µ 和 α 求二阶导数，并确认其为负值
MLE 计算示例
现在我们来看两个最大似然估计（MLE）计算的例子。
均匀分布
假设我们有一个 X 的概率分布，这意味着以下内容成立：
Pr(X=x) = 1/b-a 对于所有 a < X < b
Pr(X=x) = 0 对于所有其他 X
这里，a 和 b 是均匀分布的参数。
由于所有值的概率相同（或均匀分布），因此它被称为均匀分布。
假设数据集中有 n 个数据点，并假设这些数据点符合均匀分布。基于这些数据点，我们的目标是找到 a 和 b 的值，以定义这些数据点最有可能属于的分布。为此，我们可以使用最大似然估计方法：
我们需要找到能最大化 log(L(a,b|x[1],x[2],x[3],.....,x[n]) 的 a 和 b。
为此，我们将对 log(L(a,b|x[1],x[2],x[3],.....,x[n]) 关于 b-a 进行求导。
这给出了 -n/(b-a)，它始终小于零，表明 log(L(a,b|x[1],x[2],x[3],.....,x[n]) 是一个单调递减的函数，且其值随着 (b-a) 的增加而减少。因此，较大的 b-a 会最大化概率。
考虑到这一点，我们得到 b = max(X), a = min(X)。
泊松分布
泊松分布在本章的前面部分已有解释。简而言之，泊松分布是一个具有无限大样本数的二项分布，这么大以至于二项分布的离散性转变为了泊松分布。泊松分布也处理事件发生的概率。但不同于考虑每次试验中事件发生的概率，我们更多的是从时间间隔的角度来思考，问自己在这个时间间隔内，事件发生的次数是多少。参数也从每次试验成功的概率转变为在给定时间间隔内的成功次数。
这是总结：
-
二项分布：给定每次试验的成功概率，在给定次数的试验中取得一定次数成功的概率。
-
泊松：给定到达率或成功率的情况下，在特定时间间隔内取得一定次数成功的概率——也就是给定时间间隔内的成功平均次数。
泊松概率分布可以用以下公式表示：
这里，λ是到达率或成功率。
这个表达式给出了在给定时间间隔内观察到x次成功的概率（即到达率定义的相同时间间隔）。
我们关注的是给定一组数据集，估计符合泊松分布的λ的最大似然估计：
泊松分布的最大似然估计（MLE）计算的数学公式
注意，取对数可以在代数计算上简化计算。然而，这也引入了一些数值上的挑战——例如，确保似然值永远不为 0，因为对数不能对 0 求值。如果对数似然值无限小，数值方法也会导致无效值。
MLE 发现到达率的估计值等于数据集的均值——也就是说，过去给定时间间隔内观测到的到达次数。前面的计算可以使用 Python 中的 NumPy 和其他支持包来完成。
为了在 Python 中执行此计算，我们需要进行几个步骤：
- 编写一个函数，计算每个点的泊松概率：
import numpy as np
import math as mh
np.seterr(divide='ignore', invalid='ignore') #ignore division by zero and invalid numbers
def poissonpdf(x,lbd):
val = (np.power(lbd,x)*np.exp(-lbd))/(mh.factorial(x))
return val
- 编写一个函数，根据给定的到达率计算数据的对数似然：
def loglikelihood(data,lbd):
lkhd=1
for i in range(len(data)):
lkhd=lkhd*poissonpdf(data[i],lbd)
if lkhd!=0:
val=np.log(lkhd)
else:
val=0
return val
- 编写一个函数，计算到达率λ的对数似然的导数：
def diffllhd(data,lbd):
diff = -len(data) + sum(data)/lbd
return diff
- 生成 100 个数据点的测试数据——单位时间内的到达次数为 3 到 12 之间的随机数：
data=[randint(3, 12) for p in range(100)]
- 计算不同到达率（1 到 9）下的对数似然，并绘制图形以找到最大化的到达率：
y=[loglikelihood(data,i) for i in range(1,10)]
y=[num for num in y if num ]
x=[i for i in range(1,10) if loglikelihood(data,i)]
plt.plot(x,y)
plt.axvline(x=6,color='k')
plt.title('Log-Likelihoods for different lambdas')
plt.xlabel('Log Likelihood')
plt.ylabel('Lambda')
由此我们得到以下图形，显示当到达率为 6/单位时间时，测试数据上对数似然的最大值：
在不同的λ值（即到达率）下的对数似然值
- 使用牛顿-拉夫森方法求解对数似然的全局最大值：
def newtonRaphson(data,lbd):
h = loglikelihood(data,lbd) / diffllhd(data,lbd)
while abs(h) >= 0.0001:
if diffllhd!=0:
h = loglikelihood(data,lbd) / diffllhd(data,lbd)
# x(i+1) = x(i) - f(x) / f'(x)
lbd = lbd - h
else:
lbd=lbd
return lbd
注意：函数定义中的lbd参数是开始搜索的初始值。
牛顿-拉夫森方法是一种常用的计算方法，用于求解复杂方程的根。它是一个迭代过程，通过不断地求出自变量的不同值，直到因变量达到 0。更多信息可以参考www.math.ubc.ca/~anstee/math104/newtonmethod.pdf。
结果受初始参数值的影响很大，搜索的方向可能会因起始值的不同而大相径庭，因此在使用时要小心。
最大似然估计（MLE）概念可以扩展到执行基于分布的回归。假设我们假设到达率是一个或多个参数的函数，那么 lambda 将由这些参数的函数定义：
在这种情况下，到达率的计算方法如下：
-
在对数似然计算中使用之前方程中得到的到达率值。
-
找到对* w [0]、 w [1]和 w *[2]的对数似然的偏导数。
-
将所有偏导数等于 0，找到* w [0]、 w [1]和 w *[2]的最优值。
-
根据这些参数找到到达率的最优值。
要使用 MLE 计算，请执行以下步骤：
-
从样本参数中找到总体参数，如均值、标准差、到达率和密度。
-
对于简单线性回归无法有效拟合的数据，可以使用基于分布的回归模型，例如之前讨论的基于参数的到达率示例，或逻辑回归权重。
它在拟合回归模型中的应用使其与 OLS、梯度下降、Adam 优化、RMSprop 等优化方法处于同一类别。
参考文献
如果想更深入了解我们涉及的其他贝叶斯统计主题，请查阅以下参考文献：
-
《黑客的概率编程与贝叶斯方法》：
github.com/CamDavidsonPilon/Probabilistic-Programming-and-Bayesian-Methods-for-Hackers -
《贝叶斯数据分析》，第三版，Andrew Gelman：
www.amazon.com/Bayesian-Analysis-Chapman-Statistical-Science/dp/1439840954 -
《贝叶斯选择》，Christian P Robert（这本书更偏理论）：
www.springer.com/us/book/9780387952314
总结
本章中，我们对过去几年统计学和数据分析中最热门的趋势之一——贝叶斯统计推断方法进行了快速浏览。我们覆盖了许多内容。
我们讨论了贝叶斯统计方法的核心内容，并讨论了贝叶斯观点的多种吸引人的原因，比如它重视事实而非信仰。我们解释了关键的统计分布，并展示了如何使用各种统计包生成并在matplotlib中绘制它们。
我们在没有过度简化的情况下处理了一个相当复杂的话题，并展示了如何使用PyMC包和蒙特卡罗模拟方法，展示贝叶斯统计的强大功能，来制定模型、执行趋势分析并对真实世界数据集（Facebook 用户帖子）进行推断。最大似然估计的概念也被介绍并通过多个示例进行了说明。它是一种用于估计分布参数并将概率分布拟合到给定数据集的流行方法。
在下一章，我们将讨论如何使用 pandas 解决现实生活中的数据案例研究。
第十三章：使用 pandas 的案例研究
到目前为止，我们已经涵盖了 pandas 的广泛功能。接下来，我们将尝试在一些案例研究中实现这些功能。这些案例研究将使我们全面了解每个功能的使用，并帮助我们确定处理 DataFrame 时的关键点。此外，案例研究的逐步方法有助于加深我们对 pandas 函数的理解。本章提供了实际示例和代码片段，确保在最后，你能够理解 pandas 解决 DataFrame 问题的方法。
我们将涵盖以下案例研究：
-
从头到尾的探索性数据分析
-
使用 Python 进行网页抓取
-
数据验证
从头到尾的探索性数据分析
探索性数据分析是指理解数据特征的关键过程——如异常值、包含最相关信息的列，并通过统计和图形表示确定变量之间的关系。
让我们考虑以下 DataFrame，进行探索性数据分析：
df = pd.read_csv("data.csv")
df
以下截图展示了在 Jupyter Notebook 中加载的 DataFrame：
在 Jupyter Notebook 中加载的 DataFrame
数据概述
上述的 DataFrame 是一个汽车维修公司的客户数据。他们基本上按周期为客户提供服务。DataFrame 中的每一行对应一个独特的客户。因此，这是客户级别的数据。以下是从数据中获得的一个观察结果：
DataFrame 的形状
我们可以观察到数据包含 27,002 条记录和 26 个特征。
在开始对任何数据进行探索性数据分析之前，建议尽可能多地了解数据——包括列名及其相应的数据类型，是否包含空值（如果有，多少空值），等等。以下截图展示了通过 pandas 的info函数获得的一些基本信息：
DataFrame 的基本信息
使用info()函数，我们可以看到数据仅包含浮动和整数值。此外，没有任何列包含空值。
pandas 中的describe()函数用于获取所有数值列的各种汇总统计信息。该函数返回所有数值列的计数、均值、标准差、最小值、最大值和四分位数。以下表格展示了通过describe函数获取的数据描述：
描述数据
特征选择
如果你有一个包含多个变量的数据集，检查各列之间相关性的一个好方法是通过将相关性矩阵可视化为热图。我们可以识别并去除那些高度相关的变量，从而简化我们的分析。可视化可以通过 Python 中的seaborn库实现：
以下将是输出结果：
数据框的相关性热图
我们可以在之前的热图中观察到以下几点：
-
soldBy和days_old之间存在高度负相关 -
age_median和income_median之间存在正相关
同样地，我们可以推导出不同变量集之间的相关性。因此，基于相关性结果，我们可以通过仅选择重要特征来最小化独立特征的数量。
特征提取
除了选择有用的特征外，我们还需要从现有变量中提取显著的变量。这种方法被称为特征提取。在当前示例中，已经从现有变量中提取了一个名为new_tenure的新特征。该变量告诉我们客户在公司待了多长时间：
data['new_tenure']=data['active_status']*data['days_old']+(1-data['active_status'])*data['days_active_tenure']
以下数据框展示了新提取的变量：
含有新提取变量的数据框
数据汇总
如前所述，所呈现的数据是客户级别的数据。对汇总数据进行分析会更加可行且容易，在这种情况下，汇总数据是按区域划分的。首先，我们需要了解客户在每个区域的分布情况。因此，我们将使用groupby函数来查找每个邮政编码中的客户数量。以下代码展示了代码片段及其输出：
data.groupby('zip')['zip'].count().nlargest(10)
以下是输出结果：
基于邮政编码的汇总数据
这将给出前 10 个拥有最多客户的邮政编码。
因此，我们可以通过聚合将客户级数据转换为邮政编码级数据。在对值进行分组后，我们还必须确保去除 NA。可以使用以下代码对整个数据框进行聚合：
data_mod=data.groupby('zip')
data_clean=pd.DataFrame()
for name,data_group in data_mod:
data_group1=data_group.fillna(method='ffill')
data_clean=pd.concat([data_clean,data_group1],axis=0)
data_clean.dropna(axis=0, how='any')
以下截图是去除 NA 后的汇总数据框：
去除 NA 后的汇总数据框
data_clean将成为我们样本数据框的清理版本，该版本将传递给模型进行进一步分析。
使用 Python 进行网页抓取
网页抓取涉及从网站提取大量数据，形式可以是结构化的或非结构化的。例如，网站可能已经有一些数据以 HTML 表格元素或 CSV 文件的形式存在。这是网站上结构化数据的一个例子。但是，在大多数情况下，所需的信息会分散在网页的内容中。网页抓取有助于收集这些数据并将其存储为结构化的形式。有多种方式可以抓取网站，如在线服务、API，或者编写自己的代码。
以下是关于网页抓取的一些重要说明：
-
阅读网站的条款和条件，了解如何合法使用数据。大多数网站禁止将数据用于商业目的。
-
确保不要过快下载数据，因为这样可能会导致网站崩溃。你也有可能被网站封锁。
使用 pandas 进行网页抓取
Python 提供了不同的库来进行抓取：
-
pandas
-
BeautifulSoup
-
Scrapy
在本节中，我们将看到如何利用 pandas 和 BeautifulSoup 的强大功能进行数据抓取。首先，pandas 足以从网站上提取结构化数据，而不需要 BeautifulSoup 的帮助。在前面的章节中，我们学习了如何从不同格式（.csv、.xlsx和.xls）加载数据到 Python 中。类似于这些，pandas 有一个专门用于从 HTML 文件加载表格数据的函数。要读取 HTML 文件，pandas 的 DataFrame 会查找一个标签。这个标签称为<td> </td>标签，用于定义 HTML 中的表格。
pandas 使用read_html()来读取 HTML 文档。这个函数将 URL 中的所有结构化数据加载到 Python 环境中。因此，每当你传递一个 HTML 文件给 pandas 并希望它输出一个漂亮的 DataFrame 时，确保 HTML 页面中有一个表格。
我们可以尝试在一个示例网址上使用这个功能（www.bseindia.com/static/members/TFEquity.aspx）：
示例网页
上述网页包含多个表格。使用 pandas，我们可以提取所有表格，并将其存储在一个列表中：
包含多个 DataFrame 的列表
在下图中，正在提取网页中的第二个表格：
网页和 pandas DataFrame 的对比
清洗后，提取的 DataFrame 完全复制了网站上的内容：
清洗后的 DataFrame
通过适当的索引，所有来自网页的表格都可以通过read_html函数提取。
使用 BeautifulSoup 进行网页抓取
BeautifulSoup 是一个 Python 库（www.crummy.com/software/BeautifulSoup/），用于从 HTML 和 XML 文件中提取数据。它提供了导航、访问、搜索和修改网页 HTML 内容的方式。了解 HTML 的基础知识对于成功抓取网页内容非常重要。为了解析内容，首先我们需要做的是确定在哪里可以找到我们想要下载的文件的链接，这些文件位于 HTML 标签的多层级中。简而言之，网页上有大量代码，而我们要做的就是找到包含数据的相关代码段。
在网站上，右键点击并选择检查。这将允许你查看网站背后的原始代码。点击检查后，你应该能看到以下控制台弹出：
浏览器的检查菜单
注意我们所提到的表格被包裹在一个叫做 table 的标签中。每一行都位于 <tr> 标签之间。同样，每个单元格都位于 <td> 标签之间。理解这些基本差异可以让数据提取变得更容易。
我们首先导入以下库：
导入库
接下来，我们使用 requests 库请求 URL。如果访问成功，你应该能看到以下输出：
从网站获得成功的响应
然后，我们使用BeautifulSoup解析html，以便能够处理更整洁、嵌套的BeautifulSoup数据结构。通过一些 HTML 标签的知识，解析后的内容可以使用for循环和 pandas DataFrame 轻松转换为 DataFrame。使用 BeautifulSoup 的最大优势在于，它甚至可以从非结构化的来源中提取数据，并通过支持的库将其转化为表格，而 pandas 的read_html函数只能处理结构化的数据来源。因此，根据需求，我们使用了BeautifulSoup：
使用 BeautifulSoup 提取的 DataFrame
数据验证
数据验证是检查数据质量的过程，确保数据既正确又适用于分析。它使用称为验证规则的例程来检查输入模型的数据的真实性。在大数据时代，计算机和其他技术形式生成大量信息，这些信息推动着数据产生的数量，如果这些数据缺乏质量，那么使用它们会显得不专业，这也突出了数据验证的重要性。
在这个案例研究中，我们将考虑两个 DataFrame：
-
来自平面文件的测试 DataFrame
-
来自 MongoDB 的验证 DataFrame
在测试 DataFrame 上执行验证例程，同时将其对应的数据框作为参考。
数据概览
这里考虑的数据集是 学习管理系统（LMS）数据的一部分。它们展示了与学生注册、跟踪、报告以及教育课程的交付相关的信息。我们将从平面文件加载测试 DataFrame：
从平面文件加载测试 DataFrame
pymongo 库用于将 MongoDB 连接到 Python。通常，MongoDB 会监听端口 27017：
从 Python 连接 MongoDB
我们可以在以下截图中看到连接参数。由于数据库安装在本地，我们通过 localhost 进行连接。加载的数据库名称是 lms_db：
从 MongoDB 读取数据
结构化数据库与非结构化数据库
由于 MongoDB 属于非结构化数据库类别，因此其使用的术语与结构化数据库（如 MySQL 和 PostgreSQL）大不相同。下表展示了各种 SQL 术语和概念以及相应的 MongoDB 术语和概念：
| SQL 术语/概念 | MongoDB 术语/概念 |
|---|---|
| 数据库 | 数据库 (https://docs.mongodb.com/manual/reference/glossary/#term-database) |
| 表 | 集合 (docs.mongodb.com/manual/reference/glossary/#term-collection) |
| 行 | 文档或 BSON 文档 |
| 列 | 字段 (docs.mongodb.com/manual/reference/glossary/#term-field) |
| 索引 | 索引 (docs.mongodb.com/manual/reference/glossary/#term-index) |
| 表连接 | $lookup，嵌入式文档 (docs.mongodb.com/manual/reference/operator/aggregation/lookup/#pipe._S_lookup) |
| 主键 | 主键 (docs.mongodb.com/manual/reference/glossary/#term-primary-key) |
| 将任何唯一的列或列组合指定为主键。 | 在 MongoDB 中，主键会自动设置为_id字段。(docs.mongodb.com/manual/reference/glossary/#term-id) |
| 聚合（例如 group by） | 聚合管道 |
| 事务 | 事务 (docs.mongodb.com/manual/core/transactions/) |
SQL 与 MongoDB 术语对比视图
验证数据类型
数据类型是变量的一种属性，Python 使用它来理解如何存储和处理数据。例如，程序需要理解存储 5 和 10 的变量是数字型的，以便能够将它们相加得到 15；或者理解存储 cat 和 hat 的变量是字符串类型，以便它们可以连接（加在一起）得到 cathat。因此，它成为任何 pandas DataFrame 的初步和基本属性。
可以使用用户定义的比较函数来验证测试 DataFrame 的数据类型：
验证测试 DataFrame 的数据类型
File1 和 File2 分别对应测试数据集和验证数据集。从输出中可以明显看出，测试 DataFrame 的所有数据类型与验证 DataFrame 的数据类型匹配。如果存在不匹配，输出将显示不一致的列数。
验证维度
DataFrame 是一种二维数据结构，其中数据以表格形式呈现，类似于关系型数据库表格，按行和列排列。检查测试集和验证集是否匹配的基本方法之一是比较行数和列数。如果 DataFrame 的形状不匹配，那么测试 DataFrame 与验证 DataFrame 之间的差异就显而易见了。以下是一个截图，展示了如何验证维度：
验证维度
验证单个条目
一旦前两个测试用例满足要求，扫描单个条目以查找虚假数据就变得非常重要。前面图示中的验证过程描述了从数据采集过程中获得的值与真实值之间的差异。随着数据量的增加，验证条目变得越来越困难。通过高效使用 pandas，可以减轻这种效果。在以下示例中，使用循环（暴力法）和 pandas 索引扫描了单个条目。
使用 pandas 索引
以下截图展示了如何使用 pandas 验证单元格：
使用 pandas 索引验证单元格
使用循环
以下截图展示了如何通过使用循环验证单元格：
使用循环验证单元格
当我们使用 pandas 索引时，结果非常令人鼓舞。验证一个包含 200,000 行和 15 列的 DataFrame 只用了 0.53 秒，而使用循环完成相同的验证流程则花费了超过 7 分钟。因此，始终建议利用 pandas 的强大功能，避免使用迭代编程。
总结
pandas 对许多辅助数据活动非常有用，例如探索性数据分析、验证两个数据源之间数据的有效性（如数据类型或计数），以及构建和塑造从其他来源获取的数据，比如抓取网站或数据库。在这一章中，我们处理了这些主题的一些案例研究。数据科学家每天都会进行这些活动，本章应该能让你大致了解在真实数据集上执行这些活动的体验。
在下一章，我们将讨论 pandas 库的架构和代码结构。这将帮助我们全面了解该库的功能，并使我们能够更好地进行故障排除。
第十四章：pandas 库架构
在本章中，我们将探讨 pandas 用户可以使用的各种库。本章旨在成为一本简短的指南，帮助用户浏览和了解 pandas 提供的各种模块和库。它详细描述了库代码的组织结构，并简要介绍了各种模块。这对于希望了解 pandas 内部工作原理的用户，以及希望为代码库做出贡献的用户，将是非常有价值的。我们还将简要演示如何使用 Python 扩展来提高性能。以下是本章将讨论的各个主题：
-
pandas 库层次结构简介
-
pandas 模块和文件的描述
-
使用 Python 扩展提高性能
理解 pandas 文件层次结构
一般来说，在安装时，pandas 会作为 Python 模块安装在第三方 Python 模块的标准位置。在下表中，您可以看到 Unix/macOS 和 Windows 平台的标准安装位置：
| 平台 | 标准安装位置 | 示例 |
|---|---|---|
| Unix/macOS | prefix/lib/pythonX.Y/site-packages | /usr/local/lib/python2.7/site-packages |
| Windows | prefix\Lib\site-packages | C:\Python27\Lib\site-packages |
如果通过 Anaconda 安装 Python，那么 pandas 模块可以在 Anaconda 目录下找到，类似的文件路径为：Anaconda3\pkgs\pandas-0.23.4-py37h830ac7b_0\Lib\site-packages\pandas。
现在我们已经了解了第三方 Python 模块部分，接下来将了解文件层次结构。已安装的 Pandas 库包含八种类型的文件。这些文件遵循特定的层次结构，下面描述了这一结构：
-
pandas/core：该文件包含用于基本数据结构（如 Series/DataFrame）及相关功能的文件。 -
pandas/src：该文件包含用于实现基本算法的 Cython 和 C 代码。 -
pandas/io：该文件包含处理不同文件格式的输入/输出工具，如平面文件、Excel、HDF5 和 SQL。 -
pandas/tools：该文件包含辅助数据算法、合并和连接例程、拼接、透视表等。该模块主要用于数据操作。 -
pandas/sparse：该文件包含稀疏版本的数据结构，如系列、DataFrame、Panels 等。 -
pandas/stats：该文件包含线性回归、面板回归、移动窗口回归以及其他一些统计函数。它应该被 statsmodels 中的功能所替代。 -
pandas/util：该文件包含实用工具以及开发和测试工具。 -
pandas/rpy：该文件包含 RPy2 接口，用于连接 R，从而扩展数据分析操作的范围。
更多信息，请参见：pandas.pydata.org/developers.html。
pandas 模块和文件的描述
在本节中，我们简要描述了构成 pandas 库的各个子模块和文件。
pandas/core
该模块包含了 pandas 的核心子模块，具体讨论如下：
-
api.py：该模块导入了一些关键模块和警告信息，以供后续使用，例如索引、groupby和重塑函数。 -
apply.py：该模块包含一些类，帮助将函数应用到 DataFrame 或 series 上。 -
arrays：该模块隔离了 pandas 对numpy的依赖——即所有直接使用numpy的操作。array子模块中的base.py处理所有面向数组的操作，例如ndarray值、形状和ndim，而categorical.py子模块专门处理分类值。 -
base.py：该模块定义了诸如StringMixin和PandasObject等基础类，PandasObject是多个 pandas 对象的基类，例如Period、PandasSQLTable、sparse.array.SparseArray/SparseList、internals.Block、internals.BlockManager、generic.NDFrame、groupby.GroupBy、base.FrozenList、base.FrozenNDArray、io.sql.PandasSQL、io.sql.PandasSQLTable、tseries.period.Period、FrozenList、FrozenNDArray: IndexOpsMixin和DatetimeIndexOpsMixin。 -
common.py：该模块定义了处理数据结构的常用工具方法。例如，isnull对象用于检测缺失值。 -
config.py：该模块用于处理整个包的可配置对象。它定义了以下类：OptionError、DictWrapper、CallableDynamicDoc、option_context和config_init。 -
datetools.py：这是一个处理 Python 中日期的函数集合。它还利用了 pandas 中tseries模块的一些函数。 -
frame.py：该模块定义了 pandas 的 DataFrame 类及其各种方法。DataFrame 继承自 NDFrame（见下文）。它从pandas-core模块下的多个子模块借用函数，来定义 DataFrame 的功能性操作。 -
generic.py：该模块定义了通用的 NDFrame 基类，这是 pandas 中 DataFrame、series 和 panel 类的基类。NDFrame 继承自PandasObject，后者在base.py中定义。NDFrame 可以视为 pandas DataFrame 的 N 维版本。有关更多信息，请访问nullege.com/codes/search/pandas.core.generic.NDFrame。 -
categorical.py：该模块定义了categorical，它是从PandasObject派生的一个类，用于表示类似于 R/S-plus 的分类变量（稍后我们会进一步扩展这一点）。 -
groupby.py：该模块定义了多种类，用于实现groupby功能：-
Splitter 类：包括
DataSplitter、ArraySplitter、SeriesSplitter、FrameSplitter和NDFrameSplitter。 -
Grouper/grouping 类：包括
Grouper、GroupBy、BaseGrouper、BinGrouper、Grouping、SeriesGroupBy和NDFrameGroupBy。
-
-
ops.py：此文件定义了一个内部 API，用于对PandasObjects执行算术运算。它定义了为对象添加算术方法的函数。它定义了一个_create_methods元方法，用于通过算术、比较和布尔方法构造器创建其他方法。add_methods方法接受一个新方法列表，将它们添加到现有方法列表中，并将它们绑定到适当的类。add_special_arithmetic_methods、add_flex_arithmetic_methods、call _create_methods和add_methods用于向类中添加算术方法。
它定义了_TimeOp类，这是一个用于日期时间相关算术运算的封装类。它包含了对 series、DataFrame 和 panel 函数进行算术、比较和布尔操作的封装函数：_arith_method_SERIES(..)、_comp_method_SERIES(..)、_bool_method_SERIES(..)、_flex_method_SERIES(..)、_arith_method_FRAME(..)、_comp_method_FRAME(..)、_flex_comp_method_FRAME(..)、_arith_method_PANEL(..)和_comp_method_PANEL(..)。
-
index.py：此文件定义了索引类及其相关功能。Index 是所有 pandas 对象（如 series、DataFrame 和 panel）用来存储轴标签的工具。它下面是一个不可变的数组，提供一个有序的集合，可以进行切片操作。 -
indexing.py：此模块包含一系列函数和类，使得多重索引操作更加简便。 -
missing.py：此文件定义了诸如掩蔽和插值等技术，用于处理缺失数据。 -
internals.py：此文件定义了多个对象类，具体如下所示：-
Block：这是一个同质类型的 N 维numpy.ndarray对象，具有额外的 pandas 功能——例如，它使用__slots__来限制对象的属性为ndim、values和_mgr_locs。它作为其他 Block 子类的基类。 -
NumericBlock：这是一个用于处理数值类型区块的基类。 -
FloatOrComplexBlock：这是FloatBlock和ComplexBlock的基类，继承自`NumericBlock`。 -
ComplexBlock：这是处理复数类型区块对象的类。 -
FloatBlock：这是处理浮动类型区块对象的类。 -
IntBlock：这是处理整数类型区块对象的类。 -
TimeDeltaBlock、BoolBlock和DatetimeBlock：这些是处理时间差、布尔值和日期时间的区块类。 -
ObjectBlock：这是处理用户定义对象的区块类。 -
SparseBlock：这是处理相同类型稀疏数组的类。 -
BlockManager：这是管理一组区块对象的类，它不是公开的 API 类。 -
SingleBlockManager：这是管理单个区块的类。 -
JoinUnit：这是一个区块对象的实用类。
-
-
nanops.py：这个子模块包含一组用于专门处理 NaN 值的类和功能。 -
ops.py：该文件定义了 pandas 对象的算术运算。它不是公开 API。 -
panel.py、panel4d.py和panelnd.py：这些提供了 pandas 面板对象的功能。 -
resample.py：该文件定义了用于时间间隔分组和聚合的自定义groupby类。 -
series.py：该文件定义了 pandas Series 类及其从 NDFrame 和IndexOpsMixin继承的各种方法，以适应一维数据结构和一维时间序列数据。 -
sorting.py：该文件定义了排序所需的所有工具。 -
sparse.py：该文件定义了处理稀疏数据结构的导入。稀疏数据结构通过省略匹配 NaN 或缺失值的数据点来进行压缩。有关此的更多信息，请访问pandas.pydata.org/pandas-docs/stable/sparse.html。 -
strings.py：这些函数用于处理字符串操作，如str_replace、str_contains和str_cat。 -
window.py：该模块帮助对数据结构进行窗口处理并计算滚动窗口中的聚合值。
以下图示概述了 Pandas 核心的结构：
现在，让我们继续下一个子模块。
pandas/io
该模块包含多个数据 I/O 模块，具体如下：
-
api.py：该文件定义了数据 I/O API 的各种导入。 -
common.py：该文件定义了 I/O API 的通用功能。 -
clipboards.py：该文件包含跨平台剪贴板方法，支持通过键盘启用复制和粘贴功能。pandas I/O API 包含如pandas.read_clipboard()和pandas.to_clipboard(..)等函数。 -
date_converters.py：该文件定义了日期转换函数。 -
excel.py：该模块用于解析和转换 Excel 数据。它定义了ExcelFile和ExcelWriter类。 -
feather_format.py：该模块读取和写入 Feather 格式的数据。 -
gbq.py：这是用于 Google BigQuery 的模块。 -
html.py：这是用于处理 HTML I/O 的模块。 -
json.py：这是用于处理 pandas 中 JSON I/O 的模块。它定义了Writer、SeriesWriter、FrameWriter、Parser、SeriesParser和FrameParser类。 -
msgpack：该模块读取和写入msgpack格式的数据。 -
packer.py：该文件是一个msgpack序列化程序，支持将 pandas 数据结构读写到磁盘。 -
parquet.py：该模块读取和写入 Parquet 格式的数据。 -
parsers.py：这是定义用于解析和处理文件以创建 pandas DataFrame 的各种函数和类的模块。以下列出的三个read_*函数都有多种可配置选项用于读取。详细信息，请参见bit.ly/1EKDYbP：-
read_csv(..)：该函数定义了pandas.read_csv()，用于将 CSV 文件的内容读入 DataFrame。 -
read_table(..)：这个方法用于将制表符分隔的表文件读取到 DataFrame 中。 -
read_fwf(..)：这个方法用于将固定宽度格式的文件读取到 DataFrame 中。 -
TextFileReader：这是用于读取文本文件的类。 -
ParserBase：这是解析器对象的基类。 -
CParserWrapper、PythonParser：这些分别是 C 和 Python 的解析器。它们都继承自ParserBase。 -
FixedWidthReader：这是用于读取固定宽度数据的类。固定宽度数据文件包含在文件中特定位置的字段。 -
FixedWithFieldParser：这是用于解析固定宽度字段的类，该类继承自PythonParser。
-
-
pickle.py：该模块提供了方法来 pickle（序列化）pandas 对象，方法如下：-
to_pickle(..)：该方法用于将对象序列化到文件。 -
read_pickle(..)：这个方法用于从文件中读取序列化的对象到 pandas 对象。仅应使用受信任的源。
-
-
pytables.py：这是用于 PyTables 模块的接口，用于将 pandas 数据结构读写到磁盘上的文件。 -
sql.py：这是一个类和函数的集合，旨在从关系型数据库中检索数据，并尽量做到与数据库无关。以下是这些类和函数：-
PandasSQL：这是用于将 pandas 与 SQL 接口的基类。它提供了虚拟的read_sql和to_sql方法，这些方法必须由子类实现。 -
PandasSQLAlchemy：这是PandasSQL的子类，能够使用 SQLAlchemy 在 DataFrame 和 SQL 数据库之间进行转换。 -
PandasSQLTable：这是将 pandas 表（DataFrame）映射到 SQL 表的类。 -
pandasSQL_builder(..)：根据提供的参数返回正确的 PandasSQL 子类。 -
PandasSQLTableLegacy：这是PandasSQLTable的遗留支持版本。 -
PandasSQLLegacy：这是PandasSQLTable的遗留支持版本。 -
get_schema(..)：这个方法用于获取给定数据框架的 SQL 数据库表架构。 -
read_sql_table(..)：这个方法用于将 SQL 数据库表读取到 DataFrame 中。 -
read_sql_query(..)：这个方法用于将 SQL 查询读取到 DataFrame 中。 -
read_sql(..)：这个方法用于将 SQL 查询/表读取到 DataFrame 中。
-
-
stata.py：这个模块包含用于将 Stata 文件处理为 pandas DataFrame 的工具。 -
sas：此模块包含子模块，用于从 SAS 输出中读取数据。 -
S3.py：这个模块提供与 S3 存储桶的远程连接功能。
pandas/tools
该模块的详细信息如下：
-
plotting.py：这是用于绘图模块的包装器，最近版本中已被弃用。 -
merge.py：这个模块提供用于合并序列、DataFrame 和面板对象的函数，如merge(..)和concat(..)，并在最近的版本中已被弃用。
pandas/util
该pandas/util是提供实用功能的模块，模块的详细信息如下：
-
testing.py：该文件提供了用于测试的assertion、debug、unit test和其他类/函数。它包含许多特殊的断言函数，使得检查系列、DataFrame 或面板对象是否相等变得更加容易。一些这些函数包括assert_equal(..)、assert_series_equal(..)、assert_frame_equal(..)和assert_panelnd_equal(..)。pandas.util.testing模块对 pandas 代码库的贡献者尤其有用，它定义了一个util.TestCase类，还为潜在的代码库贡献者提供了处理区域设置、控制台调试、文件清理、比较器等的工具。 -
doctools.py：该子模块包含TablePlotter类，用于为 DataFrame 定义布局。 -
validators.py：该子模块帮助验证传递给函数的参数。例如，它帮助评估参数的长度、默认值和参数值。 -
print_versions.py：该文件定义了get_sys_info()函数，它返回系统信息字典；以及show_versions(..)函数，它显示可用 Python 库的版本。 -
misc.py：该文件定义了一些杂项工具。 -
decorators.py：该文件定义了一些装饰器函数和类。 -
替换和附加类是对函数文档字符串进行替换和附加的装饰器。有关 Python 装饰器的更多信息，请访问
www.artima.com/weblogs/viewpost.jsp?thread=240808。 -
test_decorators.py：该子模块提供了用于测试对象的装饰器。
pandas/tests
这个pandas/tests是提供 pandas 中各种对象测试的模块。具体的库文件名基本上是自解释的，这里我不再进一步详细说明；而是邀请读者自行探索。
pandas/compat
与兼容性相关的功能如下所述：
-
chainmap.py和chainmap_impl.py：这些文件提供了一个ChainMap类，可以将多个dict或映射组合在一起，产生一个可以更新的单一视图。 -
pickle_compat.py：该文件提供了在 0.12 版本之前的 pandas 对象序列化功能。
pandas/computation
这个pandas/computation是提供计算功能的模块，具体内容如下：
-
expressions.py：该文件通过numexpr提供快速的表达式计算。numexpr函数用于加速某些数值操作。它使用多个核心以及智能分块和缓存加速。它定义了evaluate(..)和where(..)方法。该模块在最新版本的 pandas 中已被弃用，替代用法将通过pandas.get_option实现。 -
有关
numexpr的更多信息，请访问code.google.com/p/numexpr/。有关此模块的使用，请访问pandas.pydata.org/pandas-docs/version/0.15.0/computation.html。
pandas/plotting
pandas/plotting 是处理所有 pandas 绘图功能的模块：
-
compat.py：此模块检查版本兼容性。 -
converter.py：此模块有助于处理用于绘图的日期时间值。它可以执行自动缩放时间序列轴和格式化日期时间轴刻度等功能。 -
core.py：此文件定义了一些帮助创建图形的类，如条形图、散点图、六边形箱形图和箱线图。 -
misc.py：此文件提供了一组绘图函数，可以将系列或 DataFrame 作为参数。此模块包含以下子模块，用于执行各种杂项任务，如绘制散点矩阵和 Andrews 曲线：-
scatter_matrix(..)：此函数绘制散点图矩阵。 -
andrews_curves(..)：此函数将多变量数据绘制为曲线，这些曲线使用样本作为傅里叶级数的系数。 -
parallel_coordinates(..)：这是一种绘图技术，可以帮助你观察数据中的聚类并直观估算统计量。 -
lag_plot(..)：此函数用于检查数据集或时间序列是否为随机的。 -
autocorrelation_plot(..)：此函数用于检查时间序列中的随机性。 -
bootstrap_plot(..)：此图用于以可视化的方式确定统计量（如均值或中位数）的不确定性。 -
radviz(..)：此图用于可视化多变量数据。
-
-
style.py：此文件提供了一组绘图样式选项。 -
timeseries.py：此文件定义了时间序列绘图的辅助类。 -
tools.py：此文件包含一些辅助函数，用于从 DataFrame 和系列创建表格布局。
pandas/tseries
本节内容涉及 pandas/tseries 模块，它赋予 pandas 处理时间序列数据的功能：
-
api.py：这是一个方便的导入集合。 -
converter.py：此文件定义了一组用于格式化和转换的类。 -
datetime：导入 pandas 后，它会通过register()函数将一组单位转换器注册到matplotlib。具体方法如下：
In [1]: import matplotlib.units as munits
In [2]: munits.registry
Out[2]: {}
In [3]: import pandas
In [4]: munits.registry
Out[4]:
{pandas.tslib.Timestamp: <pandas.tseries.converter.DatetimeConverter instance at 0x7fbbc4db17e8>,
pandas.tseries.period.Period: <pandas.tseries.converter.PeriodConverter instance at 0x7fbbc4dc25f0>,
datetime.date: <pandas.tseries.converter.DatetimeConverter instance at 0x7fbbc4dc2fc8>,
datetime.datetime: <pandas.tseries.converter.DatetimeConverter instance at 0x7fbbc4dc2a70>,
datetime.time: <pandas.tseries.converter.TimeConverter instance at 0x7fbbc4d61e18>}
-
Converter：此类包括TimeConverter、PeriodConverter和DateTimeConverter。 -
Formatters：此类包括TimeFormatter、PandasAutoDateFormatter和TimeSeries_DateFormatter。 -
Locators：此类包括PandasAutoDateLocator、MilliSecondLocator和TimeSeries_DateLocator。
Formatter 和 Locator 类用于处理 matplotlib 绘图中的刻度。
-
frequencies.py：此文件定义了指定时间序列对象频率的代码——如每日、每周、每季度、每月、每年等。此子模块依赖于 pandas/core 模块的dtypes子模块。 -
holiday.py：该文件定义了处理假期的函数和类——如Holiday、AbstractHolidayCalendar和USFederalHolidayCalendar等类。 -
offsets.py：该文件定义了各种类，包括处理与时间相关的偏移类。这些类的具体说明如下：-
DateOffset：这是一个接口，供提供时间段功能的类使用，如Week、WeekOfMonth、LastWeekOfMonth、QuarterOffset、YearOffset、Easter、FY5253和FY5253Quarter。 -
BusinessMixin：这是一个混合类，用于商业对象提供与时间相关的功能。它被BusinessDay类继承。BusinessDay子类继承自BusinessMixin和SingleConstructorOffset，并提供在工作日上的偏移。 -
MonthOffset：这是提供月时间段功能的类的接口，如MonthEnd、MonthBegin、BusinessMonthEnd和BusinessMonthBegin。 -
MonthEnd和MonthBegin：这两个提供了在月末或月初的日期偏移，偏移量为一个月。 -
BusinessMonthEnd和BusinessMonthBegin：这两个提供了在商业日历的月末或月初的日期偏移，偏移量为一个月。 -
YearOffset：这个偏移量是由提供年周期功能的类继承的——如YearEnd、YearBegin、BYearEnd和BYearBegin。 -
YearEnd和YearBegin：这两个提供了在年末或年初的日期偏移，偏移量为一年。 -
BYearEnd和BYearBegin：这两个提供了在商业日历的结束或开始时的日期偏移，偏移量为一年。 -
Week：它提供了一个星期的偏移。 -
WeekDay：它提供了从星期几（例如，Tue）到一周中的任何一天（例如，=2）的映射。 -
WeekOfMonth和LastWeekOfMonth：这描述了每个月中某周的日期。 -
QuarterOffset：这是提供季度周期功能的类的基类——如QuarterEnd、QuarterBegin、BQuarterEnd和BQuarterBegin。 -
QuarterEnd、QuarterBegin、BQuarterEnd和BQuarterBegin：这些类与Year*类似，不同之处在于其周期是季度，而不是年份。 -
FY5253和FY5253Quarter：这些类分别描述了 52 周和 53 周的财政年度，也称为 4-4-5 日历。 -
Easter：这是DateOffset类，用于复活节假期。 -
Tick：这是时间单位类的基类，例如Day、Hour、Minute、Second、Milli、Micro和Nano。
-
-
plotting.py：该文件从pandas-plotting模块导入tsplot(..)子模块。
接下来，我们将看到如何使用 Python 扩展来提高 Python 代码的性能。
使用 Python 扩展提高性能
Python 和 pandas 用户的一个抱怨是，虽然语言和模块的易用性和表达性非常强大，但也伴随着显著的缺点——性能。这种情况尤其在数值计算时显现。
根据编程基准标准，Python 在许多算法或数据结构操作上通常比编译语言（如 C/C++）慢。例如，在一个模拟实验中，Python3 的运行速度比最快的 C++ 实现的 n 体模拟计算慢了 104 倍。
那么，我们如何解决这个合法却令人烦恼的问题呢？我们可以在保持 Python 中我们喜爱的东西——清晰性和生产力——的同时，减缓 Python 的运行速度。我们可以通过编写代码中的性能敏感部分（例如，数字处理、算法等）为 C/C++ 代码，并通过编写 Python 扩展模块让 Python 调用这些代码来实现。更多详细信息，请访问 docs.python.org/2/extending/extending.html。
Python 扩展模块使我们能够从 Python 调用用户定义的 C/C++ 代码或库函数，从而提高代码性能，并享受 Python 使用的便利。
为了帮助我们理解 Python 扩展模块是什么，让我们考虑在 Python 中导入一个模块时发生了什么。导入语句导入一个模块，但这究竟意味着什么呢？有三种可能性，具体如下：
-
一些 Python 扩展模块在构建解释器时与之链接。
-
导入语句会导致 Python 加载
.pyc文件到内存中。.pyc文件包含 Python 的字节码，如以下代码片段所示：
In [3]: import pandas
pandas.__file__
Out[3]: '/usr/lib/python2.7/site-packages/pandas/__init__.pyc'
- 导入语句会导致 Python 扩展模块加载到内存中。
.so（共享对象）文件包含了机器码，如以下代码片段所示：
In [4]: import math
math.__file__
Out[4]: '/usr/lib/python2.7/lib-dynload/math.so'
我们将重点讨论第三种可能性，因为这是最常见的情况。即使我们处理的是从 C 编译出来的二进制共享对象，我们仍然可以将其作为 Python 模块导入。这展示了 Python 扩展的强大——应用程序可以从 Python 机器码或机器码导入模块，接口是相同的。Cython 和 SWIG 是用 C 和 C++ 编写扩展的两种最流行的方法。在编写扩展时，我们将 C/C++ 机器码封装起来，并将其转化为表现得像纯 Python 代码的 Python 扩展模块。在这次简短的讨论中，我们将只关注 Cython，因为它是专为 Python 设计的。
Cython 是 Python 的一个超集，旨在通过允许我们调用外部编译的 C/C++ 代码，并声明变量的类型，显著提高 Python 的性能。
Cython 命令从 Cython 源文件生成优化的 C/C++ 源文件，并将此优化后的 C/C++ 源文件编译成 Python 扩展模块。它为 NumPy 提供了内置支持，并将 C 的性能与 Python 的可用性结合起来。
我们将快速演示如何使用 Cython 显著加速我们的代码。让我们定义一个简单的斐波那契函数：
In [17]: def fibonacci(n):
a,b=1,1
for i in range(n):
a,b=a+b,a
return a
In [18]: fibonacci(100)
Out[18]: 927372692193078999176L
In [19]: %timeit fibonacci(100)
100000 loops, best of 3: 18.2 µs per loop
使用 timeit 模块，我们发现每次循环需要 18.2 微秒。
现在让我们通过以下步骤，在 Cython 中重写函数，并为变量指定类型：
- 首先，我们在 iPython 中导入 Cython 魔法函数，如下所示：
In [22]: %load_ext cythonmagic
- 接下来，我们在 Cython 中重写我们的函数，并为我们的变量指定类型：
In [24]: %%cython
def cfibonacci(int n):
cdef int i, a,b
for i in range(n):
a,b=a+b,a
return a
- 让我们测试一下我们新的 Cython 函数的执行时间：
In [25]: %timeit cfibonacci(100)
1000000 loops, best of 3: 321 ns per loop
In [26]: 18.2/0.321
Out[26]: 56.69781931464174
- 我们可以看到，Cython 版本比纯 Python 版本快了 57 倍！
有关使用 Cython/SWIG 或其他选项编写 Python 扩展的更多信息，请参考以下来源：
-
Pandas 文档，标题为 提升性能，链接：
pandas.pydata.org/pandas-docs/stable/enhancingperf.html -
ScipPy 讲义，标题为 与 C 接口，链接：
scipy-lectures.github.io/advanced/interfacing_with_c/interfacing_with_c.html -
Cython 文档，链接：
docs.cython.org/index.html -
SWIG 文档，链接：
www.swig.org/Doc2.0/SWIGDocumentation.html
总结
总结本章内容，我们浏览了 pandas 库的层次结构，试图说明该库的内部构造。理解这些内容对于从 pandas 代码构建自定义模块或作为开源贡献者改进 pandas 的功能将非常有帮助。我们还讨论了通过使用 Python 扩展模块来加速代码性能的好处。
在下一章中，我们将看到 pandas 与其他数据分析工具在各种分析操作方面的比较。
第十五章：pandas 与其他工具的比较
本章重点比较 pandas 与 R（许多 pandas 功能的模型工具）、SQL 和 SAS 等工具的异同，后者与 pandas 有着显著的重叠。本章旨在为希望使用 pandas 的 R、SQL 和 SAS 用户提供指南，也为希望在 pandas 中重现其代码功能的用户提供帮助。本章重点介绍 R、SQL 和 SAS 用户可用的一些关键功能，并通过一些示例演示如何在 pandas 中实现类似功能。本章假设您已安装 R 统计软件包。如果没有，您可以从此处下载并安装：www.r-project.org/。
在本章结束时，数据分析用户应该能够很好地掌握这些工具相对于 pandas 的数据分析能力，从而在需要时能够顺利过渡到或使用 pandas。工具比较的各个因素包括：
-
数据类型及其 pandas 对应物
-
切片与选择
-
数据类型列的算术运算
-
聚合与 GroupBy
-
匹配
-
分割-应用-合并
-
数据重塑与熔化
-
因子和分类数据
与 R 的比较
R 是 pandas 设计灵感的工具。两者在语法、用法和输出方面非常相似。主要的差异出现在某些数据类型上，例如 R 中的矩阵与 pandas 中的数组，R 中的 aggregate 函数与 pandas 中的 GroupBy 操作，以及一些功能相似的函数在语法上的细微差别，如 melt 和 cut。
R 中的数据类型
R 具有五种原始类型或原子类型：
-
字符
-
数值型
-
整数
-
复数
-
逻辑/布尔值
它还具有以下更复杂的容器类型：
-
向量：这与
numpy.array相似。它只能包含相同类型的对象。 -
列表：这是一个异构容器。它在 pandas 中的对应物是一个序列（series）。
-
数据框（DataFrame）：这是一个异构的二维容器，相当于 pandas 的 DataFrame。
-
矩阵：这是一个同质的二维版本的向量。它类似于
numpy.array。
本章我们将重点介绍列表和数据框，它们在 pandas 中的对应物是：序列（series）和数据框（DataFrame）。
有关 R 数据类型的更多信息，请参考以下文档：www.statmethods.net/input/datatypes.html。
有关 NumPy 数据类型的更多信息，请参考以下文档：docs.scipy.org/doc/numpy/reference/generated/numpy.array.html 和 docs.scipy.org/doc/numpy/reference/generated/numpy.matrix.html。
R 列表
R 列表可以通过显式的列表声明创建，如下所示：
>h_lst<- list(23,'donkey',5.6,1+4i,TRUE)
>h_lst
[[1]]
[1] 23
[[2]]
[1] "donkey"
[[3]]
[1] 5.6
[[4]]
[1] 1+4i
[[5]]
[1] TRUE
>typeof(h_lst)
[1] "list"
以下代码块展示了其在 pandas 中的系列等效版本，其中包括列表的创建，随后从中创建 Series：
In [8]: h_list=[23, 'donkey', 5.6,1+4j, True]
In [9]: import pandas as pd
h_ser=pd.Series(h_list)
In [10]: h_ser
Out[10]: 0 23
1 donkey
2 5.6
3 (1+4j)
4 True
dtype: object
在 pandas 中，数组的索引从 0 开始，而在 R 中，索引从 1 开始。以下是一个示例：
In [11]: type(h_ser)
Out[11]: pandas.core.series.Series
R DataFrames
我们可以通过调用data.frame()构造函数来构造 R DataFrame，然后如下面所示展示它：
>stocks_table<- data.frame(Symbol=c('GOOG','AMZN','FB','AAPL',
'TWTR','NFLX','LINKD'),
Price=c(518.7,307.82,74.9,109.7,37.1,
334.48,219.9),
MarketCap=c(352.8,142.29,216.98,643.55,23.54,20.15,27.31))
>stocks_table
Symbol PriceMarketCap
1 GOOG 518.70 352.80
2 AMZN 307.82 142.29
3 FB 74.90 216.98
4 AAPL 109.70 643.55
5 TWTR 37.10 23.54
6 NFLX 334.48 20.15
7 LINKD 219.90 27.31
在以下代码块中，我们构造了一个 pandas DataFrame 并展示出来：
In [29]: stocks_df=pd.DataFrame({'Symbol':['GOOG','AMZN','FB','AAPL',
'TWTR','NFLX','LNKD'],
'Price':[518.7,307.82,74.9,109.7,37.1,
334.48,219.9],
'MarketCap($B)' : [352.8,142.29,216.98,643.55,
23.54,20.15,27.31]
})
stocks_df=stocks_df.reindex_axis(sorted(stocks_df.columns,reverse=True),axis=1)
stocks_df
Out[29]:
Symbol PriceMarketCap($B)
0 GOOG 518.70 352.80
1 AMZN 307.82 142.29
2 FB 74.90 216.98
3 AAPL 109.70 643.55
4 TWTR 37.10 23.54
5 NFLX 334.48 20.15
6 LNKD219.90 27.31
切片和选择
在 R 中，我们通过以下三种方式进行对象切片：
-
[: 这总是返回与原始对象相同类型的对象，并且可以用来选择多个元素。 -
[[：用于提取列表或 DataFrame 的元素，只能用于提取单个元素。返回的元素类型不一定是列表或 DataFrame。 -
$：用于通过名称提取列表或 DataFrame 的元素，类似于[[。
以下是 R 中一些切片示例及其在 pandas 中的等效形式：
比较 R 矩阵和 NumPy 数组
让我们看一下 R 中的创建和选择：
>r_mat<- matrix(2:13,4,3)
>r_mat
[,1] [,2] [,3]
[1,] 2 6 10
[2,] 3 7 11
[3,] 4 8 12
[4,] 5 9 13
要选择第一行，我们写出如下代码：
>r_mat[1,]
[1] 2 6 10
要选择第二列，我们使用如下命令：
>r_mat[,2]
[1] 6 7 8 9
现在，我们来看 NumPy 数组的创建和选择：
In [60]: a=np.array(range(2,6))
b=np.array(range(6,10))
c=np.array(range(10,14))
In [66]: np_ar=np.column_stack([a,b,c])
np_ar
Out[66]: array([[ 2, 6, 10],
[ 3, 7, 11],
[ 4, 8, 12],
[ 5, 9, 13]])
要选择第一行，我们使用如下命令：
In [79]: np_ar[0,]
Out[79]: array([ 2, 6, 10])
R 和 pandas/NumPy 的索引方式不同。
在 R 中，索引从 1 开始，而在 pandas/NumPy 中，索引从 0 开始。因此，在将 R 转换为 pandas/NumPy 时，我们需要将所有索引减去 1。
要选择第二列，我们使用如下命令：
In [81]: np_ar[:,1]
Out[81]: array([6, 7, 8, 9])
另一种方法是先转置数组，然后选择列，如下所示：
In [80]: np_ar.T[1,]
Out[80]: array([6, 7, 8, 9])
比较 R 列表和 pandas 系列
在 R 中，列表的创建和选择如下所示：
>cal_lst<- list(weekdays=1:8, mth='jan')
>cal_lst
$weekdays
[1] 1 2 3 4 5 6 7 8
$mth
[1] "jan"
>cal_lst[1]
$weekdays
[1] 1 2 3 4 5 6 7 8
>cal_lst[[1]]
[1] 1 2 3 4 5 6 7 8
>cal_lst[2]
$mth
[1] "jan"
在 pandas 中，Series 的创建和选择如下所示：
In [92]: cal_df= pd.Series({'weekdays':range(1,8), 'mth':'jan'})
In [93]: cal_df
Out[93]: mthjan
weekdays [1, 2, 3, 4, 5, 6, 7]
dtype: object
In [97]: cal_df[0]
Out[97]: 'jan'
In [95]: cal_df[1]
Out[95]: [1, 2, 3, 4, 5, 6, 7]
In [96]: cal_df[[1]]
Out[96]: weekdays [1, 2, 3, 4, 5, 6, 7]
dtype: object
在这里，我们可以看到 R 列表和 pandas 系列在[]和[[]]运算符下的不同。通过考虑第二个项目——一个字符字符串，我们可以看出区别。
在 R 中，[]运算符返回容器类型，即一个包含字符串的列表，而[[]]运算符返回原子类型，在这个例子中是字符类型，如下所示：
>typeof(cal_lst[2])
[1] "list"
>typeof(cal_lst[[2]])
[1] "character"
在 pandas 中，情况正好相反：[]返回原子类型，而[[]]返回复杂类型，即 Series，如下所示：
In [99]: type(cal_df[0])
Out[99]: str
In [101]: type(cal_df[[0]])
Out[101]: pandas.core.series.Series
在 R 和 pandas 中，都可以通过列名来指定元素。
在 R 中指定列名
在 R 中，可以通过在列名前加上$运算符来完成，如下所示：
>cal_lst$mth
[1] "jan"
> cal_lst$'mth'
[1] "jan"
在 pandas 中指定列名
在 pandas 中，我们可以像往常一样使用方括号中的列名来获取子集：
In [111]: cal_df['mth']
Out[111]: 'jan'
R 和 pandas 在嵌套元素的子集操作上有所不同。例如，要从工作日中获取第四天，我们必须在 R 中使用[[]]运算符：
>cal_lst[[1]][[4]]
[1] 4
>cal_lst[[c(1,4)]]
[1] 4
然而，在 pandas 中，我们只需使用双[]：
In [132]: cal_df[1][3]
Out[132]: 4
R DataFrame 与 pandas DataFrame
在 R DataFrame 和 pandas DataFrame 中选择数据遵循类似的脚本。以下部分解释了我们如何从两者中进行多列选择。
R 中的多列选择
在 R 中，我们通过在方括号内指定向量来选择多个列：
>stocks_table[c('Symbol','Price')]
Symbol Price
1 GOOG 518.70
2 AMZN 307.82
3 FB 74.90
4 AAPL 109.70
5 TWTR 37.10
6 NFLX 334.48
7 LINKD 219.90
>stocks_table[,c('Symbol','Price')]
Symbol Price
1 GOOG 518.70
2 AMZN 307.82
3 FB 74.90
4 AAPL 109.70
5 TWTR 37.10
6 NFLX 334.48
7 LINKD 219.90
在 pandas 中的多列选择
在 pandas 中，我们按照通常的方式使用列名进行子集选择，即列名放在方括号中：
In [140]: stocks_df[['Symbol','Price']]
Out[140]:Symbol Price
0 GOOG 518.70
1 AMZN 307.82
2 FB 74.90
3 AAPL 109.70
4 TWTR 37.10
5 NFLX 334.48
6 LNKD 219.90
In [145]: stocks_df.loc[:,['Symbol','Price']]
Out[145]: Symbol Price
0 GOOG 518.70
1 AMZN 307.82
2 FB 74.90
3 AAPL 109.70
4 TWTR 37.10
5 NFLX 334.48
6 LNKD 219.90
列上的算术操作
在 R 和 pandas 中，我们可以以类似的方式在数据列中应用算术运算。因此，我们可以对两个或多个 DataFrame 中对应位置的元素执行加法或减法等算术操作。
在这里，我们构建一个 R 中的 DataFrame，列标为 x 和 y，并将列 y 从列 x 中减去：
>norm_df<- data.frame(x=rnorm(7,0,1), y=rnorm(7,0,1))
>norm_df$x - norm_df$y
[1] -1.3870730 2.4681458 -4.6991395 0.2978311 -0.8492245 1.5851009 -1.4620324
R 中的 with 运算符也具有与算术操作相同的效果：
>with(norm_df,x-y)
[1] -1.3870730 2.4681458 -4.6991395 0.2978311 -0.8492245 1.5851009 -1.4620324
在 pandas 中，相同的算术操作可以在列上执行，相应的运算符是 eval：
In [10]: import pandas as pd
import numpy as np
df = pd.DataFrame({'x': np.random.normal(0,1,size=7), 'y': np.random.normal(0,1,size=7)})
In [11]: df.x-df.y
Out[11]: 0 -0.107313
1 0.617513
2 -1.517827
3 0.565804
4 -1.630534
5 0.101900
6 0.775186
dtype: float64
In [12]: df.eval('x-y')
Out[12]: 0 -0.107313
1 0.617513
2 -1.517827
3 0.565804
4 -1.630534
5 0.101900
6 0.775186
dtype: float64
聚合与 GroupBy
有时，我们希望将数据拆分成子集，并对每个子集应用一个函数，如平均值、最大值或最小值。在 R 中，我们可以通过 aggregate 或 tapply 函数来实现。
在这里，我们有一个数据集，包含了 2014 年欧洲冠军联赛足球赛半决赛四支球队中五位前锋的统计数据。我们将使用它来说明 R 中的聚合及其在 pandas 中的等效 GroupBy 功能。
R 中的聚合
在 R 中，聚合是通过以下命令实现的：
> goal_stats=read.csv('champ_league_stats_semifinalists.csv')
>goal_stats
Club Player Goals GamesPlayed
1 Atletico Madrid Diego Costa 8 9
2 Atletico Madrid ArdaTuran 4 9
3 Atletico Madrid RaúlGarcía 4 12
4 Atletico Madrid AdriánLópez 2 9
5 Atletico Madrid Diego Godín 2 10
6 Real Madrid Cristiano Ronaldo 17 11
7 Real Madrid Gareth Bale 6 12
8 Real Madrid Karim Benzema 5 11
9 Real Madrid Isco 3 12
10 Real Madrid Ángel Di María 3 11
11 Bayern Munich Thomas Müller 5 12
12 Bayern Munich ArjenRobben 4 10
13 Bayern Munich Mario Götze 3 11
14 Bayern Munich Bastian Schweinsteiger 3 8
15 Bayern Munich Mario Mandzukić 3 10
16 Chelsea Fernando Torres 4 9
17 Chelsea Demba Ba 3 6
18 Chelsea Samuel Eto'o 3 9
19 Chelsea Eden Hazard 2 9
20 Chelsea Ramires 2 10
我们现在计算每个前锋的每场比赛进球比率，以此来衡量他们在门前的致命性：
>goal_stats$GoalsPerGame<- goal_stats$Goals/goal_stats$GamesPlayed
>goal_stats
Club Player Goals GamesPlayedGoalsPerGame
1 Atletico Madrid Diego Costa 8 9 0.8888889
2 Atletico Madrid ArdaTuran 4 9 0.4444444
3 Atletico Madrid RaúlGarcía 4 12 0.3333333
4 Atletico Madrid AdriánLópez 2 9 0.2222222
5 Atletico Madrid Diego Godín 2 10 0.2000000
6 Real Madrid Cristiano Ronaldo 17 11 1.5454545
7 Real Madrid Gareth Bale 6 12 0.5000000
8 Real Madrid Karim Benzema 5 11 0.4545455
9 Real Madrid Isco 3 12 0.2500000
10 Real Madrid Ángel Di María 3 11 0.2727273
11 Bayern Munich Thomas Müller 5 12 0.4166667
12 Bayern Munich ArjenRobben 4 10 0.4000000
13 Bayern Munich MarioGötze 3 11 0.2727273
14 Bayern Munich Bastian Schweinsteiger 3 8 0.3750000
15 Bayern Munich MarioMandzukić 3 10 0.3000000
16 Chelsea Fernando Torres 4 9 0.4444444
17 Chelsea Demba Ba 3 6 0.5000000
18 Chelsea Samuel Eto'o 3 9 0.3333333
19 Chelsea Eden Hazard 2 9 0.2222222
20 Chelsea Ramires 2 10 0.2000000
假设我们想知道每个球队的最高进球比率。我们可以通过以下方式计算：
>aggregate(x=goal_stats[,c('GoalsPerGame')], by=list(goal_stats$Club),FUN=max)
Group.1 x
1 Atletico Madrid 0.8888889
2 Bayern Munich 0.4166667
3 Chelsea 0.5000000
4 Real Madrid 1.5454545
tapply 函数用于对由一个或多个列定义的数组或向量的子集应用函数。tapply 函数也可以如下使用：
>tapply(goal_stats$GoalsPerGame,goal_stats$Club,max)
Atletico Madrid Bayern Munich Chelsea Real Madrid
0.8888889 0.4166667 0.5000000 1.5454545
pandas 中的 GroupBy 运算符
在 pandas 中，我们可以通过使用 GroupBy 函数来实现相同的结果：
In [6]: import pandas as pd
importnumpy as np
In [7]: goal_stats_df=pd.read_csv('champ_league_stats_semifinalists.csv')
In [27]: goal_stats_df['GoalsPerGame']= goal_stats_df['Goals']/goal_stats_df['GamesPlayed']
In [27]: goal_stats_df['GoalsPerGame']= goal_stats_df['Goals']/goal_stats_df['GamesPlayed']
In [28]: goal_stats_df
Out[28]: Club Player Goals GamesPlayedGoalsPerGame
0 Atletico Madrid Diego Costa 8 9 0.888889
1 Atletico Madrid ArdaTuran 4 9 0.444444
2 Atletico Madrid RaúlGarcía 4 12 0.333333
3 Atletico Madrid AdriánLópez 2 9 0.222222
4 Atletico Madrid Diego Godín 2 10 0.200000
5 Real Madrid Cristiano Ronaldo 17 11 1.545455
6 Real Madrid Gareth Bale 6 12 0.500000
7 Real Madrid Karim Benzema 5 11 0.454545
8 Real Madrid Isco 3 12 0.250000
9 Real Madrid Ángel Di María 3 11 0.272727
10 Bayern Munich Thomas Müller 5 12 0.416667
11 Bayern Munich ArjenRobben 4 10 0.400000
12 Bayern Munich Mario Götze 3 11 0.272727
13 Bayern Munich BastianSchweinsteiger 3 8 0.375000
14 Bayern Munich MarioMandzukić 3 10 0.300000
15 Chelsea Fernando Torres 4 9 0.444444
16 Chelsea Demba Ba 3 6 0.500000
17 Chelsea Samuel Eto'o 3 9 0.333333
18 Chelsea Eden Hazard 2 9 0.222222
19 Chelsea Ramires 2 10 0.200000
In [30]: grouped = goal_stats_df.groupby('Club')
In [17]: grouped['GoalsPerGame'].aggregate(np.max)
Out[17]: Club
Atletico Madrid 0.888889
Bayern Munich 0.416667
Chelsea 0.500000
Real Madrid 1.545455
Name: GoalsPerGame, dtype: float64
In [22]: grouped['GoalsPerGame'].apply(np.max)
Out[22]: Club
Atletico Madrid 0.888889
Bayern Munich 0.416667
Chelsea 0.500000
Real Madrid 1.545455
Name: GoalsPerGame, dtype: float64
比较 R 和 pandas 中的匹配运算符
在这里，我们演示了 R 中的匹配运算符 (%in%) 和 pandas 中的匹配运算符 (isin()) 之间的等效性。在这两种情况下，都会生成一个逻辑向量（R）或系列（pandas），它表示找到匹配的位置信息。
R 中的 %in% 运算符
在这里，我们演示了如何在 R 中使用 %in% 运算符：
>stock_symbols=stocks_table$Symbol
>stock_symbols
[1] GOOG AMZN FB AAPL TWTR NFLX LINKD
Levels: AAPL AMZN FB GOOG LINKD NFLX TWTR
>stock_symbols %in% c('GOOG','NFLX')
[1] TRUE FALSE FALSE FALSE FALSE TRUE FALSE
pandas 中的 isin() 函数
下面是一个使用 pandas isin() 函数的例子：
In [11]: stock_symbols=stocks_df.Symbol
stock_symbols
Out[11]: 0 GOOG
1 AMZN
2 FB
3 AAPL
4 TWTR
5 NFLX
6 LNKD
Name: Symbol, dtype: object
In [10]: stock_symbols.isin(['GOOG','NFLX'])
Out[10]: 0 True
1 False
2 False
3 False
4 False
5 True
6 False
Name: Symbol, dtype: bool
逻辑子集选择
在 R 和 pandas 中，有多种方式可以执行逻辑子集选择。假设我们希望显示所有平均每场比赛进球数大于或等于 0.5 的球员；也就是说，平均每两场比赛至少进一球。
R 中的逻辑子集选择
这是我们在 R 中如何实现的：
- 使用逻辑切片：
>goal_stats[goal_stats$GoalsPerGame>=0.5,]
Club Player Goals GamesPlayedGoalsPerGame
1 Atletico Madrid Diego Costa 8 9 0.8888889
6 Real Madrid Cristiano Ronaldo 17 11 1.5454545
7 Real Madrid Gareth Bale 6 12 0.5000000
17 Chelsea Demba Ba 3 6 0.5000000
- 使用
subset()函数：
>subset(goal_stats,GoalsPerGame>=0.5)
Club Player Goals GamesPlayedGoalsPerGame
1 Atletico Madrid Diego Costa 8 9 0.8888889
6 Real Madrid Cristiano Ronaldo 17 11 1.5454545
7 Real Madrid Gareth Bale 6 12 0.5000000
17 Chelsea Demba Ba 3 6 0.5000000
pandas 中的逻辑子集选择
在 pandas 中，我们做类似的操作：
- 逻辑切片：
In [33]: goal_stats_df[goal_stats_df['GoalsPerGame']>=0.5]
Out[33]: Club Player Goals GamesPlayedGoalsPerGame
0 Atletico Madrid Diego Costa 8 9 0.888889
5 Real Madrid Cristiano Ronaldo 17 11 1.545455
6 Real Madrid Gareth Bale 6 12 0.500000
16 Chelsea Demba Ba 3 6 0.500000
DataFrame.query()运算符：
In [36]: goal_stats_df.query('GoalsPerGame>= 0.5')
Out[36]:
Club Player Goals GamesPlayedGoalsPerGame
0 Atletico Madrid Diego Costa 8 9 0.888889
5 Real Madrid Cristiano Ronaldo 17 11 1.545455
6 Real Madrid Gareth Bale 6 12 0.500000
16 Chelsea Demba Ba 3 6 0.500000
分割-应用-合并
R 有一个叫做plyr的库，用于拆分-应用-合并数据分析。plyr库有一个函数叫做ddply，它可以用来将一个函数应用到 DataFrame 的子集上，然后将结果合并成另一个 DataFrame。
有关ddply的更多信息，请参考以下链接：www.inside-r.org/packages/cran/plyr/docs/ddply。
为了说明，让我们考虑一个最近创建的 R 数据集的子集，包含 2013 年从纽约市出发的航班数据：cran.r-project.org/web/packages/nycflights13/index.html。
在 R 中的实现
在这里，我们在 R 中安装该包并实例化库：
>install.packages('nycflights13')
...
>library('nycflights13')
>dim(flights)
[1] 336776 16
>head(flights,3)
year month day dep_timedep_delayarr_timearr_delay carrier tailnum flight
1 2013 1 1 517 2 830 11 UA N14228 1545
2 2013 1 1 533 4 850 20 UA N24211 1714
3 2013 1 1 542 2 923 33 AA N619AA 1141
origindestair_time distance hour minute
1 EWR IAH 227 1400 5 17
2 LGA IAH 227 1416 5 33
3 JFK MIA 160 1089 5 42
> flights.data=na.omit(flights[,c('year','month','dep_delay','arr_delay','distance')])
>flights.sample<- flights.data[sample(1:nrow(flights.data),100,replace=FALSE),]
>head(flights.sample,5)
year month dep_delayarr_delay distance
155501 2013 3 2 5 184
2410 2013 1 0 4 762
64158 2013 11 -7 -27 509
221447 2013 5 -5 -12 184
281887 2013 8 -1 -10 937
ddply函数使我们能够按年份和月份总结出发延迟（均值和标准差）：
>ddply(flights.sample,.(year,month),summarize, mean_dep_delay=round(mean(dep_delay),2), s_dep_delay=round(sd(dep_delay),2))
year month mean_dep_delaysd_dep_delay
1 2013 1 -0.20 2.28
2 2013 2 23.85 61.63
3 2013 3 10.00 34.72
4 2013 4 0.88 12.56
5 2013 5 8.56 32.42
6 2013 6 58.14 145.78
7 2013 7 25.29 58.88
8 2013 8 25.86 59.38
9 2013 9 -0.38 10.25
10 2013 10 9.31 15.27
11 2013 11 -1.09 7.73
12 2013 12 0.00 8.58
让我们将flights.sample数据集保存为 CSV 文件，这样我们就可以使用这些数据向我们展示如何在 pandas 中做相同的操作：
>write.csv(flights.sample,file='nycflights13_sample.csv', quote=FALSE,row.names=FALSE)
在 pandas 中的实现
为了在 pandas 中做同样的操作，我们读取前面章节保存的 CSV 文件：
In [40]: flights_sample=pd.read_csv('nycflights13_sample.csv')
In [41]: flights_sample.head()
Out[41]: year month dep_delayarr_delay distance
0 2013 3 2 5 184
1 2013 1 0 4 762
2 2013 11 -7 -27 509
3 2013 5 -5 -12 184
4 2013 8 -1 -10 937
我们通过使用GroupBy()操作符实现了与ddply相同的效果，如下所示的代码和输出：
In [44]: pd.set_option('precision',3)
In [45]: grouped = flights_sample_df.groupby(['year','month'])
In [48]: grouped['dep_delay'].agg([np.mean, np.std])
Out[48]: mean std
year month
2013 1 -0.20 2.28
2 23.85 61.63
3 10.00 34.72
4 0.88 12.56
5 8.56 32.42
6 58.14 145.78
7 25.29 58.88
8 25.86 59.38
9 -0.38 10.25
10 9.31 15.27
11 -1.09 7.73
12 0.00 8.58
使用 melt 进行重塑
melt函数将宽格式数据转换为由唯一 ID-变量组合组成的单列数据。
R 中的 melt 函数
这里，我们展示了在 R 中使用melt()函数。它生成了长格式数据，其中每一行都是独特的变量-值组合：
>sample4=head(flights.sample,4)[c('year','month','dep_delay','arr_delay')]
> sample4
year month dep_delayarr_delay
155501 2013 3 2 5
2410 2013 1 0 4
64158 2013 11 -7 -27
221447 2013 5 -5 -12
>melt(sample4,id=c('year','month'))
year month variable value
1 2013 3 dep_delay 2
2 2013 1 dep_delay 0
3 2013 11 dep_delay -7
4 2013 5 dep_delay -5
5 2013 3 arr_delay 5
6 2013 1 arr_delay 4
7 2013 11 arr_delay -27
8 2013 5 arr_delay -12
>
更多信息请参考以下链接：www.statmethods.net/management/reshape.html。
pandas 的 melt 函数
在 pandas 中，melt函数类似：
In [55]: sample_4_df=flights_sample_df[['year','month','dep_delay', \
'arr_delay']].head(4)
In [56]: sample_4_df
Out[56]: year month dep_delayarr_delay
0 2013 3 2 5
1 2013 1 0 4
2 2013 11 -7 -27
3 2013 5 -5 -12
In [59]: pd.melt(sample_4_df,id_vars=['year','month'])
Out[59]: year month variable value
0 2013 3 dep_delay 2
1 2013 1 dep_delay 0
2 2013 11 dep_delay -7
3 2013 5 dep_delay -5
4 2013 3 arr_delay 5
5 2013 1 arr_delay 4
6 2013 11 arr_delay -27
7 2013 5 arr_delay -12
该信息的参考来源如下：pandas.pydata.org/pandas-docs/stable/reshaping.html#reshaping-by-melt。
分类数据
在 R 中，分类变量称为因子，cut()函数使我们能够将连续的数值变量分割成不同的范围，并将这些范围视为因子或分类变量，或者将分类变量划分到更大的区间中。
R 示例使用 cut()
以下代码块展示了 R 中的一个示例：
clinical.trial<- data.frame(patient = 1:1000,
age = rnorm(1000, mean = 50, sd = 5),
year.enroll = sample(paste("19", 80:99, sep = ""),
1000, replace = TRUE))
>clinical.trial<- data.frame(patient = 1:1000,
+ age = rnorm(1000, mean = 50, sd = 5),
+ year.enroll = sample(paste("19", 80:99, sep = ""),
+ 1000, replace = TRUE))
>summary(clinical.trial)
patient age year.enroll
Min. : 1.0 Min. :31.14 1995 : 61
1st Qu.: 250.8 1st Qu.:46.77 1989 : 60
Median : 500.5 Median :50.14 1985 : 57
Mean : 500.5 Mean :50.14 1988 : 57
3rd Qu.: 750.2 3rd Qu.:53.50 1990 : 56
Max. :1000.0 Max. :70.15 1991 : 55
(Other):654
>ctcut<- cut(clinical.trial$age, breaks = 5)> table(ctcut)
ctcut
(31.1,38.9] (38.9,46.7] (46.7,54.6] (54.6,62.4] (62.4,70.2]
15 232 558 186 9
pandas 解决方案
以下代码块包含了之前解释的cut()函数在 pandas 中的等效实现（仅适用于版本 0.15 及以上）：
In [79]: pd.set_option('precision',4)
clinical_trial=pd.DataFrame({'patient':range(1,1001),
'age' : np.random.normal(50,5,size=1000),
'year_enroll': [str(x) for x in np.random.choice(range(1980,2000),size=1000,replace=True)]})
In [80]: clinical_trial.describe()
Out[80]: age patient
count 1000.000 1000.000
mean 50.089 500.500
std 4.909 288.819
min 29.944 1.000
25% 46.572 250.750
50% 50.314 500.500
75% 53.320 750.250
max 63.458 1000.000
In [81]: clinical_trial.describe(include=['O'])
Out[81]: year_enroll
count 1000
unique 20
top 1992
freq 62
In [82]: clinical_trial.year_enroll.value_counts()[:6]
Out[82]: 1992 62
1985 61
1986 59
1994 59
1983 58
1991 58
dtype: int64
In [83]: ctcut=pd.cut(clinical_trial['age'], 5)
In [84]: ctcut.head()
Out[84]: 0 (43.349, 50.052]
1 (50.052, 56.755]
2 (50.052, 56.755]
3 (43.349, 50.052]
4 (50.052, 56.755]
Name: age, dtype: category
Categories (5, object): [(29.91, 36.646] < (36.646, 43.349] < (43.349, 50.052] < (50.052, 56.755] < (56.755, 63.458]]
In [85]: ctcut.value_counts().sort_index()
Out[85]: (29.91, 36.646] 3
(36.646, 43.349] 82
(43.349, 50.052] 396
(50.052, 56.755] 434
(56.755, 63.458] 85
dtype: int64
前面章节中的比较可以通过以下表格进行总结：
R 和 pandas 中数据结构与操作的比较
与 SQL 的比较
pandas 在许多方面与 SQL 类似，它用于数据选择、数据过滤、数据聚合、数据生成和数据修改。SQL 对数据库表执行的操作类似于 pandas 对 DataFrame 所执行的操作。在本节中，我们将比较 SQL 中的功能与其在 pandas 中的等效功能。
SELECT
SELECT 用于选择或子集化表中某些列的数据。假设你有一个名为 DallasData 的表/DataFrame。这些数据可以附在你的书包中，或者从书中的云盘访问。要从三个给定列中选择五行数据，你可以写出如下命令：
SQL
在 SQL 中，你可以使用以下命令：
select state_name,active_status,services_due from DallasData LIMIT 5;
pandas
在 pandas 中，你可以使用以下命令：
DallasData[['state_name','active_status','services_due']].head(5)
以下是前述命令的输出结果：
DallasData 上的选择结果
Where
Where 语句在 SQL 中用于应用过滤条件，根据特定标准筛选行。在 pandas 中，相应的操作是基于条件的逻辑子集选择。
假设我们想找出 active_status == 1 的行。这可以在这两种工具中如下进行。
SQL
在 SQL 中，你可以使用以下命令：
select * from DallasData where active_status ==1 LIMIT 5;
pandas
在 pandas 中，你可以使用以下命令：
DallasData[DallasData['active_status']==1].head(5);
以下是前述命令的输出结果：
过滤掉只有活跃客户后的 DallasData
假设我们想找出活跃的客户（active_status == 1）且已完成的服务少于九项（services_completed < 9）的行。这可以在这两种工具中如下进行。
SQL
在 SQL 中，你可以使用以下命令：
select * from DallasData where active_status ==1 AND services_completed <9 LIMIT 5;
pandas
在 pandas 中，你可以使用以下命令：
DallasData[(DallasData['active_status']==1) & (DallasData['services_completed'] <9)].head(5)
以下是前述命令的输出结果：
过滤掉活跃且已完成超过九项服务的客户后的 DallasData
假设我们想找出活跃的客户（active_status == 1）的行，但只查找这些行的客户 ID、邮政编码和卖家 ID。这可以在这两种工具中如下进行。
SQL
在 SQL 中，你可以使用以下命令：
select customerID,zip,soldBy from DallasData where active_status ==1 LIMIT 5;
pandas
在 pandas 中，你可以使用以下命令：
DallasData[DallasData['active_status']==1][['customerID','zip','soldBy']].head(5)
以下是前述命令的输出结果：
过滤掉只有活跃客户且只选择特定列后的 DallasData
group by
group by 语句用于聚合数据，并查找数值列的聚合值。执行此操作的关键字相同，但语法略有不同。我们来看看几个示例。
假设我们要查找数据集中活跃和非活跃客户的数量。这可以在这两种工具中如下进行。
SQL
在 SQL 中，你可以使用以下命令：
select active_status, count(*) as number from DallasData group by active_status;
pandas
在 pandas 中，你可以使用以下命令：
DallasData.groupby('active_status').size();
以下是前述命令的输出结果：
使用 Python 中的 groupby 统计活跃和非活跃客户的数量
不同的聚合操作可以同时应用于两个不同的列，下面的示例展示了如何操作。
SQL
在 SQL 中，你可以使用以下命令：
select active_status, sum(services_complted), mean(age_median) from DallasData group by active_status;
pandas
在 pandas 中，你可以使用以下命令：
DallasData.groupby('active_status').agg({'services_completed':np.sum,'age_median':np.mean})
按照活跃和非活跃客户分组，统计已完成的服务总数和平均客户年龄，使用 Python 中的 groupby
对多个列进行聚合或多重索引聚合也是可行的。假设我们希望按邮政编码获取活跃和非活跃客户的详细信息。
SQL
在 SQL 中，你可以使用以下命令：
select active_status, sum(services_complted), mean(days_old) from DallasData group by active_status,zip;
pandas
在 pandas 中，你可以使用以下命令：
DallasData.groupby(['active_status','zip']).agg({'services_completed':np.sum,'days_old':np.mean}).head(5)
以下是前述命令的输出：
使用 Python 中的 groupby 按客户的活跃状态和邮政编码进行多重索引分组
update
SQL 中的update语句用于根据某些条件筛选数据行，并更新或修改这些行中的某些值。在 pandas 中，没有特定的关键字或函数来执行此操作；相反，它是通过直接赋值来完成的。让我们看几个例子。
假设已经确定数据管理员在数据收集过程中犯了错误。由于这个错误，实际上年龄为 45 的数据点被随机赋值为大于 35 的数值。为了解决这个问题，我们将把所有这些行（age>35）更新为 45。
SQL
在 SQL 中，你可以使用以下命令：
update DallasData set age_median=45 where age_median>35
pandas
在 pandas 中，你可以使用以下命令：
DallasData[DallasData['age_median']>35]=45
以下两个截图展示了执行更新操作前后的数据：
更新所有年龄大于 35 岁为 45 岁之前和之后的数据
delete
SQL 中的delete语句用于根据某些条件从数据库表中删除数据行。在 pandas 中，我们并不删除行，而是选择取消选中它们。让我们来看一些示例。
假设我们想查看那些在系统中至少存在 500 天的客户（days_old>500）。
SQL
在 SQL 中，你可以使用以下命令：
delete DallasData where days_old<500
pandas
在 pandas 中，你可以使用以下命令：
DallasData1 = DallasData[DallasData['days_old']>500]
运行以下命令以检查是否执行了预期的操作。
DallasData1[DallasData1['days_old']<400]
如果删除操作执行正确，这应返回 0 行数据。
JOIN
join语句用于合并数据库中的不同表，并提取分散在多个表中的重要信息。在 pandas 中，merge 操作符完成了相同的工作。唯一的区别是语法略有不同。
让我们创建两个数据集，用来演示 SQL 和 pandas 中的不同连接及其语法：
df1 = pd.DataFrame({'key': ['IN', 'SA', 'SL', 'NZ'],'Result':['W','L','L','W']})
df2 = pd.DataFrame({'key': ['IN', 'SA', 'SA', 'WI'],'Score':[200,325,178,391]})
以下是前述命令的输出：
两个示例数据集
假设我们想在这两者之间进行内连接。可以按照前述的方法在这两个工具中实现。
SQL
在 SQL 中，你可以使用以下命令：
SELECT * FROM df1 INNER JOIN df2 ON df1.key = df2.key;
pandas
在 pandas 中，你可以使用以下命令：
pd.merge(df1,df2,on='key')
以下是前面命令的输出：
两个 DataFrame 的内连接输出
正如内连接所预期的，只有在两个表中都存在的键值才会出现在合并后的数据集中。
假设我们要在两者之间实现左连接。这可以通过以下两种工具来完成：
SQL
在 SQL 中，你可以使用以下命令：
SELECT * FROM df1 LEFT JOIN df2 ON df1.key = df2.key;
pandas
在 pandas 中，你可以使用以下命令：
pd.merge(df1,df2,on='key',how='left')
以下是前面命令的输出：
两个表的左连接输出
正如左连接所预期的，它会检索左表（此处为 df1）中存在的所有唯一键值及其在右表中的相应值。对于左表中的键值，如果在右表中没有找到匹配项，则返回 NaN。
假设我们要在两者之间实现右连接。这可以通过以下两种工具来完成：
SQL
在 SQL 中，你可以使用以下命令：
SELECT * FROM df1 RIGHT JOIN df2 ON df1.key = df2.key;
pandas
在 pandas 中，你可以使用以下命令：
pd.merge(df1,df2,on='key',how='right')
以下是前面命令的输出：
两个表的右连接输出
正如右连接所预期的，它会检索右表（此处为 df2）中存在的所有唯一键值及其在左表中的相应值。对于右表中的键值，如果在左表中没有找到匹配项，则返回 NaN。
与 SAS 的比较
SAS 是过去的分析利器。在 R 和 Python 这两个开源运动的代表工具取代它之前，SAS 是分析解决方案的市场领导者，曾居于头号地位。然而，尽管其成本过高，许多企业仍然依赖它处理所有的分析需求。
在本节中，我们将所有比较都以表格的形式展示。SAS 和 pandas 的对应关系总结在下表中：
| Pandas | SAS |
|---|---|
| DataFrame | dataset |
| column | variable |
| row | observation |
| groupby | BY-group |
| NaN | . |
现在，让我们看看如何在 pandas 和 SAS 中执行基本的数据操作：
| 任务 | Pandas | SAS |
|---|---|---|
| 创建数据集 | pd.DataFrame({'odds': [1, 3, 5, 7, 9], 'evens': [2, 4, 6, 8, 10]}) | data df; input x y; datalines; 1 2 3 4 5 6 7 8 9 10; run; |
| 读取数据集 | pd.read_csv('DallasData.csv') | proc import datafile='DallasData.csv' dbms=csv out=tips replace; getnames=yes; run; |
| 导出数据集 | DallasData.to_csv('dallas.csv') | proc export data=DallasData outfile='dallas.csv' dbms=csv; run; |
| 列操作 | DallasData['days_old_year'] = DallasData['days_old']/365 | data DallasData;`` set DallasData;`` days_old_year = days_old / 365;``run; |
| 过滤 | DallasData[DallasData['days_old']>800].head() | data tips;`` set DallasData;`` if days_old > 800;``run; |
| If-else | DallasData['income_class'] = np.where(DallasData['income_average'] < 40000, 'low', 'high') | data DallasData;`` set dallas;`` format income_average $5.;`` if days_old < 40000 then bucket = 'low';`` else bucket = 'high';``run; |
| 列选择 | DallasData[['zip','customerID','days_old','services_due']].head() | data dallas;`` set DallasData;`` keep zip CustomerID days_old services_due;``run; |
| 排序 | dallas = DallasData.sort_values(['days_old','services_completed']) | proc sort data=DallasData;`` by days_old services_completed;``run; |
| 字符串长度 | DallasData['state_name'].str.len().head() | data _null_;``set DallasData;``put(LENGTHN(state_name));``put(LENGTHC(state_name));``run; |
| 分组聚合 | dallas_grouped = DallasData.groupby(['zip', 'customerID'])['days_old', 'services_completed'].sum() | proc summary data=DallasData nway;`` class zip customerID;`` var days_old services_completed;`` output out=dallas_summed sum=;``run; |
| 联接 | df1 = pd.DataFrame({'key': ['A', 'B', 'C', 'D'], 'value': np.random.randn(4)})``df2 = pd.DataFrame({'key': ['B', 'D', 'D', 'E'],'value': np.random.randn(4)})``inner_join = df1.merge(df2, on=['key'], how='inner')``left_join = df1.merge(df2, on=['key'], how='left')``right_join = df1.merge(df2, on=['key'], how='right') | proc sort data=df1;`` by key;``run;``proc sort data=df2;`` by key;``run;``data left_join inner_join right_join outer_join;`` merge df1(in=a) df2(in=b);`` if a and b then output inner_join;`` if a then output left_join;`` if b then output right_join;`` if a or b then output outer_join;``run; |
总结
在本章中，我们尝试将 R 和 SQL 中的关键功能与 pandas 中的对应功能进行比较，以实现以下目标：
-
帮助 R、SQL 和 SAS 用户，他们可能希望在 pandas 中复制相同的功能
-
为了帮助任何阅读 R、SQL 和 SAS 代码的用户，他们可能希望将代码重写为 pandas 代码
在下一章中，我们将通过简要介绍 scikit-learn 库来进行机器学习，并演示 pandas 如何融入该框架，从而总结本书。本章的参考文档可以在此找到：pandas.pydata.org/pandas-docs/stable/comparison_with_r.html。
