5.4 实战：假设检验与空间聚类

5.4 实战：假设检验与空间聚类

本节摘要：一个实战小节，两半任务。前半是统计推断的规范流程：两组实验数据先做描述统计与正态性检验，再根据方差齐性选择独立样本 t 检验或曼惠特尼 U 检验，最后计算效应量并给出结论。后半是空间分析：用 KDTree 对城市坐标做最近邻查询，用凸包观察城市群布局。全程可运行的代码，串起前两节的知识。

学习目标

阅读完本节，你应当能够：

1. 独立走完"描述统计 → 正态性检验 → 方差齐性 → 选检验 → 效应量 → 结论"的完整假设检验流程；
2. 在 t 检验与曼惠特尼 U 检验之间做出有依据的选择，并解释理由；
3. 用 KDTree 完成城市点集的最近邻与半径查询，并解读结果；
4. 用 ConvexHull 分析点集的整体边界，把空间结论写进报告。

一、问题与直觉

一个实验，两组数据，领导要结论。看起来只是"算个 p 值"，但直接甩一个 t 检验出去，等于把前提假设当空气。真实的科研流程里，检验方法的选择本身就是要论证的一环：数据正态吗？方差齐吗？不满足前提的检验，p 值再小也没人信。

这一节我们不只跑代码，还跑流程。任务有两个：第一，比较两种肥料的作物产量差异，给出"有没有差异、差异多大、结论怎么说"的完整回答；第二，给一份城市坐标，找出每个城市最近的邻居，并判断这些城市在空间上是否聚成几簇。第一个任务用统计，第二个任务用空间算法——正好把 5.1 与 5.2 的知识各用一遍。

二、实战一：假设检验完整流程

2.1 数据与描述统计

先造两组模拟数据：对照组产量均值一百，实验组一百零八，标准差都是十二，各三十个样本。真实项目里这一步是读数据、清洗、画箱线图，这里直接生成：

import numpy as np
from scipy import stats

rng = np.random.default_rng(42)
ctrl = rng.normal(100, 12, size=30)
treat = rng.normal(108, 12, size=30)

print(stats.describe(ctrl))
print(stats.describe(treat))
​

描述统计给出两组各自的均值、方差、偏度、峰度。先看数字再往下走：均值差约八个点，方差量级接近，偏度都不大——初步印象是"有点差异，前提可能满足"。注意，这是印象，不是结论。

2.2 正态性检验

t 检验的前提是两组数据各自近似正态。用 shapiro 检验逐组验证：

print(stats.shapiro(ctrl))
print(stats.shapiro(treat))
​

两个 p 值都远大于 0.05，不能拒绝"数据正态"的零假设。这里有一个容易犯的认知错误：p 值大于 0.05 不等于"确认正态"，只能说没有证据说它不正态。样本量小时，shapiro 的检验力有限；但三十个样本配合偏度峰度一起看，够用了。

2.3 方差齐性检验

t 检验还要求两组方差接近。用 levene 检验：

print(stats.levene(ctrl, treat))
​

p 值远大于 0.05，方差齐性假设成立。到这一步，两条前提都过了，可以放心用标准 t 检验；如果 levene 检验显著，就要改用 equal_var=False 的 Welch 校正版本。

2.4 选择检验并执行

前提满足，走参数检验；不满足，走非参数检验。执行：

res = stats.ttest_ind(ctrl, treat, equal_var=True)
print(res)                       # 统计量与 p 值
​

p 值小于 0.05，结论：两组产量差异在 0.05 显著性水平下统计显著。作为对照，我们再看如果数据是偏态的，曼惠特尼 U 检验怎么用：

skewed_a = rng.exponential(10, size=30)
skewed_b = rng.exponential(13, size=30)
print(stats.shapiro(skewed_a).pvalue)      # 非正态
print(stats.mannwhitneyu(skewed_a, skewed_b).pvalue)
​

曼惠特尼 U 检验不假设分布形态，只比较两组数据的秩，是非正态场景下的默认选择。

2.5 效应量：显著之外的另一半答案

p 值回答"差异是不是抽样波动"，效应量回答"差异有多大"。t 检验的常用效应量是 Cohen d：两组均值差除以合并标准差。

n1, n2 = len(ctrl), len(treat)
sp = np.sqrt(((n1 - 1) * ctrl.var(ddof=1) + (n2 - 1) * treat.var(ddof=1)) / (n1 + n2 - 2))
cohen_d = (treat.mean() - ctrl.mean()) / sp
print(round(cohen_d, 3))
​

按经验标准，d 约 0.2 是小效应，0.5 是中效应，0.8 是大效应。报告里写"p 小于 0.05，效应量 d 约 0.65"，比只写 p 值有信息量得多——样本量翻倍，p 值会变小，但效应量基本不变，这才是差异的真实大小。

图 5-3 假设检验决策流程

2.6 结论怎么写

报告不能只有一行"p 小于 0.05"。完整结论要包含四件事：用了什么检验及前提验证结果、p 值与效应量、置信信息、实际意义。表 5-2 是这套流程的输出清单：

2.7 数据是怎么来的：模拟数据与真实流程

这一节从头到尾用的都是模拟数据，这不是偷懒，而是刻意为之：只有数据生成过程完全已知，才能验证整个检验流程对不对。真值摆在那里（均值差八、标准差十二），流程跑完如果结论对不上，说明哪一步写错了。真实项目里，这一步对应的是数据读取、缺失值处理、离群值筛查和箱线图检查——流程一样，只是数据来源不同。另一个要点是随机种子：default_rng(42) 固定了随机序列，代码在任何机器上重跑结果一致。改掉种子，数字会变，但结论方向不会变，这正是抽样波动的含义：单次结果不可复现的实验，结论没人敢信。

2.8 置信区间：给结论画一个边界

p 值给出"是否显著"的二元判断，置信区间给出"差异可能有多大"的范围，两者互补。两组各三十个样本、均值差八、合并标准差十二的设计，95% 置信区间大约从 1.8 到 14.2：

n = 30
diff, sp = 8.0, 12.0
se = sp * np.sqrt(1 / n + 1 / n)
lo = diff - stats.t.ppf(0.975, 2 * n - 2) * se
hi = diff + stats.t.ppf(0.975, 2 * n - 2) * se
print(lo, hi)      # 约 1.8 与 14.2
​

区间不含零，与 p 值结论一致；区间宽度则揭示了不确定性——下限 1.8 意味着真实差异可能只是小幅改善，上限 14.2 意味着可能接近设计值的两倍。报告里写"均值差约八，95% 置信区间 1.8 到 14.2"，比单写一个 p 值信息量大得多。注意样本量减半时区间会明显变宽，这就是为什么功效不足的实验常常给不出有意义的结论。

2.9 功效：这个实验设计够不够

按设计参数算一下功效：效应量零点六七（八除以十二）、每组三十个样本、显著性水平零点零五，功效约七成。也就是说，即使差异真实存在，这个实验仍有近三成概率检不出。想达到九成功效，每组样本量要提到五十左右。如果项目只够三十个样本，那结论的措辞就要克制："未检出显著差异"与"不存在差异"是两回事。功效分析应该在实验设计阶段做，而不是等结果出来再补——跑完再算功效，属于事后诸葛亮，参考价值有限。

三、实战二：空间聚类与邻近分析

3.1 城市坐标的最近邻

给八个城市的经纬度，找出每个城市最近的城市。坐标数量少，暴力也能算，但我们用 KDTree 走标准流程，因为真实项目里可能是几百个站点：

from scipy.spatial import KDTree, ConvexHull

cities = np.array([
    [116.4, 39.9],   # 北京
    [121.5, 31.2],   # 上海
    [113.3, 23.1],   # 广州
    [114.1, 22.5],   # 深圳
    [104.1, 30.7],   # 成都
    [108.9, 34.3],   # 西安
    [117.2, 39.1],   # 天津
    [118.8, 32.1],   # 南京
])
tree = KDTree(cities)
for i, city in enumerate(cities):
    dist, idx = tree.query(city, k=2)   # k=1 是自己，k=2 是最近邻
    print(i, dist[1], idx[1])
​

k=2 是因为查询点本身一定排第一，第二个才是真正的最近邻。结果里北京与天津互为最近邻，上海与南京互为最近邻，广州与深圳互为最近邻——三对"双子城"当场现形。这就是最近邻分析的价值：不需要任何先验，从数据里自己长出结构。

再试试半径查询：以北京为中心，半径 1.5 度内有哪些城市？

print(tree.query_ball_point(cities[0], r=1.5))
​

返回的列表里除了北京自己，还有天津。半径查询在气象上叫"邻近站点检索"：第 6 章的综合案例里，找与某站点相关性最强的邻近站点，用的就是同一套机制。

3.2 凸包看布局

最近邻回答局部关系，凸包回答整体边界：

hull = ConvexHull(cities)
print(hull.vertices)
​

八个城市里，位于边缘的五个（北京、天津、上海、广州、深圳）构成凸包，成都与西安落在内部。这个结果直观反映了中国东部城市的空间分布：沿海一线是边界，内陆城市在包内。聚类观察到这里可以画图验证，但结论已经清楚：城市群在空间上聚成东、南、西南几簇，KDTree 的邻近查询给出了定量证据。

3.3 把结果画出来：邻近关系的可视化解读

最近邻列表是数字，画出来才直观。把每个城市和它的最近邻连一条线：

import matplotlib.pyplot as plt

fig, ax = plt.subplots()
ax.plot(cities[:, 0], cities[:, 1], 'o')
for i in range(len(cities)):
    dist, idx = tree.query(cities[i], k=2)
    ax.plot([cities[i, 0], cities[idx[1], 0]],
            [cities[i, 1], cities[idx[1], 1]], color='gray', lw=0.8)
​

图上一眼能看到四对互连关系：沿海的北京天津、上海南京、广州深圳三对双子城，内陆还有成都与西安互为最近邻。最近邻分析的价值就在这里：不给任何先验，结构自己长出来。再画凸包，会有一个反直觉的细节：西安与成都也在凸包顶点里，南京反而在内部。别急着怀疑算法——这是经纬度直接当平面坐标的几何失真，纬线在高纬度收缩，最西端的点被拉到了边界上。要做地理上正确的边界分析，先把经纬度投影到等距或等积投影坐标系再算，结论才靠得住。

3.4 两条数据流的汇合

两个任务表面无关，骨架相同：先造数据、再检查前提、然后分析、最后解读。整条管线画出来是两路并进、最后汇成一份报告：

统计那边输出"差异是否显著、有多大"，空间那边输出"谁离谁近、边界在哪"，合起来就是一份完整的数据分析报告。第 6 章的气象案例会沿用同一套骨架：数据清洗、前提检查、算法分析、结论解读，四步一个不少。

四、工程实践要点

两个任务走下来，几条经验值得固化。第一，检验流程的顺序不能跳：先正态、再齐性、后主检验，每一步的输出都是下一步选择的依据。第二，报告里 p 值与效应量缺一不可，只报 p 值等于只报一半结论。第三，空间分析里 k=2 的"自己排第一"细节容易踩，批量查询前先拿单点验证一遍。第四，经纬度直接当平面坐标用，在几百公里的尺度上误差可接受，但做高精度分析要投影到平面坐标系。

⚠️ 常见坑：shapiro 检验对样本量敏感，样本太大时微小偏差也会显著，样本太小时检验力不足——结合偏度峰度和直方图一起判断，别只信一个 p 值。KDTree 查询时忘记 k 里包含查询点自己，最近邻下标取错，是空间分析里最常见的低级错误。
💡 关键直觉：统计检验与空间分析看似两个世界，本质都是"从数据里找结构"。前者找的是差异结构，后者找的是位置结构；共同的纪律是先检查前提，再下结论。

五、常见失败与修正

两个任务跑完，把最容易翻车的地方列成一张表，对照着自查：

最后一条最常见也最隐蔽：检验数量一多，纯随机就能撞出显著结果，二十个检验里出现一个 p 小于 0.05 完全正常。先想清楚要回答几个问题，再决定做几个检验，而不是跑完再挑好看的报。

⚠️ 常见坑：模拟数据的结论别直接搬进真实报告——真数据的生成机制未知，抽样误差、测量误差、混杂因素都在，结论的表述要比模拟案例谨慎得多。

核心回顾

• 要点一：假设检验的规范顺序是描述统计、正态性、方差齐性、主检验、效应量，一步不能跳。
• 要点二：正态且方差齐用 t 检验，正态但方差不齐用 Welch 校正，非正态用曼惠特尼 U 检验。
• 要点三：效应量 Cohen d 不随样本量膨胀，报告差异大小时比 p 值更诚实。
• 要点四：KDTree 查询时 k=1 返回查询点自己，取最近邻要用 k=2 或过滤下标。
• 要点五：半径查询 query_ball_point 是邻近站点检索的标准工具，第 6 章会复用。
• 要点六：凸包给出点集的整体边界，内部点与边界点一望便知。
• 要点七：完整结论包含检验名称、前提验证、p 值、效应量与实际意义，缺一不可。

到此，统计、空间、特殊函数三块拼图就位。第 6 章把它们和前面所有章节装进同一条气象数据分析管线，那是全书最后一次大串联。