统计分析

统计分析

概述

进行假设检验（t检验、方差分析、卡方检验）、回归、相关性和贝叶斯分析，并系统地检查假设、计算效应量并撰写APA风格的报告。目标是产出一份评审者挑不出毛病的分析：正确的检验方法、经过验证的假设、如实报告的效应量，以及完整的撰写。

何时使用此技能

在以下情况下使用此技能：

• 进行统计假设检验（t检验、方差分析、卡方检验、非参数检验）
• 执行回归或相关性分析
• 运行贝叶斯统计分析
• 检查统计假设和诊断
• 计算效应量并进行功效分析
• 以APA格式报告统计结果
• 分析研究的实验或观察数据

安装

使用 uv 安装本技能所用的库。生产环境中请固定版本号；探索阶段不固定版本也可以。

# 核心频率派工具栈（Python 3.10+；推荐3.12+以使用最新的SciPy/ArviZ）
uv pip install "pingouin>=0.6" "scipy>=1.11" "statsmodels>=0.14.6" pandas matplotlib seaborn

# 贝叶斯建模（PyMC 5 + ArviZ）
uv pip install "pymc>=5.0" "arviz>=1.0"

兼容性说明（已针对 pingouin 0.6.1、statsmodels 0.14.6、arviz 1.2 验证，2026年）：

• Pingouin 0.6.0 重命名了输出列名以去除特殊字符：p_val、cohen_d、CI95、p_unc（在0.5.x中此前分别为p-val、cohen-d、CI95%、p-unc）。下文示例使用当前列名；如果仍在使用0.5.x，请改用带连字符的旧写法。
• statsmodels + SciPy：使用 statsmodels>=0.14.6 搭配 scipy>=1.11，以避免在SciPy 1.16+上出现_lazywhere导入错误。
• ArviZ 1.x：az.summary() 现在默认输出89%区间（eti89列），且区间宽度参数改为 ci_prob（而不是hdi_prob）。要报告常规的95%可信区间，需显式传入 az.summary(trace, ci_prob=0.95)。
• 单侧贝叶斯因子已从Pingouin中移除：pg.ttest(..., alternative='greater') 会静默丢弃 BF10 列，pg.bayesfactor_ttest 在单侧备择假设下会报错。对于单侧贝叶斯检验，请直接使用PyMC（计算方向性假设的后验概率）或JASP/R的BayesFactor包。

有关特定模型的API（OLS、GLM、ARIMA），请参见 statsmodels 技能。有关PyMC工作流，请参见 pymc 技能。

分析工作流

每一次严谨的分析都遵循相同的路径。跳过步骤正是分析最终被撤稿的原因，因此请按顺序逐一完成，并说明在每一步做了什么。

1. 先明确问题，再动数据。 说明假设、结局变量和预测变量，以及设计类型（独立还是配对、组数）。现在就确定计划采用的检验方法——在看过结果之后再选择检验方法，即使是无意为之，也是p值操纵（p-hacking）。
2. 检查数据。 按组统计：n、均值、标准差、中位数、缺失值。在进行任何检验之前先绘制原始数据（直方图或箱线图）。组间样本量不均、缺失数据、地板/天花板效应和异常值都会改变适用的检验方法——应向用户明确指出这些问题，而不是悄悄地绕过它们。
3. 选择检验方法，可使用下面的快速参考，或参考references/test_selection_guide.md了解基础方法之外的设计（计数、时间-事件、信度、析因设计）。
4. 用`scripts/assumption_checks.py`检查假设。 如果某个假设不成立，改用补救性检验方法（见下表），并同时报告原计划和所做的更改。
5. 运行检验，并始终同时计算效应量——p值只能说明效应是否存在；效应量才能说明是否值得关注。
6. 撰写报告，使用下面的APA模板，包括描述性统计、精确的统计量、带置信区间的效应量，以及所执行的假设检查。

如果用户只需要其中一个步骤（例如"我需要多少参与者？"），可以直接跳到相应部分——但仍需确认该计算所依赖的设计假设。

检验选择指南

快速参考：选择正确的检验

有关计数数据、生存分析、信度、析因设计等全面指导，请使用references/test_selection_guide.md。快速参考：

比较两组：

• 独立、连续、正态 → 独立t检验
• 独立、连续、非正态 → 曼-惠特尼U检验
• 配对、连续、正态 → 配对t检验
• 配对、连续、非正态 → 威尔科克森符号秩检验
• 二分类结果 → 卡方检验或Fisher精确检验

比较3+组：

• 独立、连续、正态 → 单因素方差分析
• 独立、连续、非正态 → 克鲁斯卡尔-沃利斯检验
• 配对、连续、正态 → 重复测量方差分析
• 配对、连续、非正态 → 弗里德曼检验

关系：

• 两个连续变量 → 皮尔逊（正态）或斯皮尔曼相关（非正态）
• 连续结果与预测变量 → 线性回归
• 二分类结果与预测变量 → 逻辑回归

贝叶斯替代方法：

所有检验都有贝叶斯版本，可提供关于假设的直接概率陈述、量化证据的贝叶斯因子，以及支持零假设的能力。参见references/bayesian_statistics.md。

假设检查

始终在解释检验结果之前检查假设，并报告检查过程——评审者会关注这些内容。

使用内置的scripts/assumption_checks.py模块。从技能目录（skills/statistical-analysis/）运行Python，或将scripts/添加到sys.path：

from assumption_checks import comprehensive_assumption_check

# 异常值 + 正态性（按组）+ 方差齐性检验，附带图表
results = comprehensive_assumption_check(
    data=df,
    value_col='score',
    group_col='group',  # 可选：用于组比较
    alpha=0.05
)

对于有针对性的检查，可导入单独的函数：

from assumption_checks import (
    check_normality,                # 夏皮罗-威尔克检验 + Q-Q图 + 直方图
    check_normality_per_group,
    check_homogeneity_of_variance,  # 列文检验 + 箱线图
    check_linearity,                # 简单回归的散点图 + 残差图
    check_regression_diagnostics,   # 完整的OLS诊断（见下文回归部分）
    detect_outliers                 # IQR或z分数方法
)

result = check_normality(data=df['score'], name='Test Score', alpha=0.05, plot=True)
print(result['interpretation'])
print(result['recommendation'])

假设被违反时的处理

正态性被违反：

• 轻微违反 + 每组n > 30 → 继续使用参数检验（稳健）
• 中度违反 → 使用非参数替代方法
• 严重违反 → 转换数据或使用非参数检验

方差齐性被违反：

• 对于t检验 → 使用Welch t检验（pg.ttest在correction='auto'下会自动应用）
• 对于方差分析 → 使用Welch方差分析（pg.welch_anova）或Brown-Forsythe检验
• 对于回归 → 使用稳健标准误或加权最小二乘法

线性被违反（回归）：

• 添加多项式项、转换变量，或使用非线性模型/GAM

正式检验会随着n的增大而变得过于敏感：当n ≥ 100时，应更重视Q-Q图而不是夏皮罗-威尔克检验的p值。有关全面指导，请参见references/assumptions_and_diagnostics.md。

运行统计检验

主要使用的库：

• pingouin：用户友好的检验方法，默认返回效应量——标准检验优先使用它
• scipy.stats：核心统计检验
• statsmodels：回归、诊断、功效分析
• pymc + arviz：贝叶斯建模和诊断

带完整报告的T检验

import pingouin as pg

# correction='auto' 在方差不齐时会自动应用Welch校正
result = pg.ttest(group_a, group_b, correction='auto')

# Pingouin >= 0.6 的列名
t_stat = result['T'].values[0]
df = result['dof'].values[0]
p_value = result['p_val'].values[0]
cohens_d = result['cohen_d'].values[0]
ci_lower, ci_upper = result['CI95'].values[0]  # 均值差的置信区间

print(f"t({df:.0f}) = {t_stat:.2f}, p = {p_value:.3f}, d = {cohens_d:.2f}")

带事后检验的方差分析

import pingouin as pg

aov = pg.anova(dv='score', between='group', data=df, detailed=True)
print(aov)

# 效应量：部分η²
eta_p2 = aov['np2'].values[0]

# 如果显著，进行事后检验（Tukey HSD控制族错误率）
if aov['p_unc'].values[0] < 0.05:
    posthoc = pg.pairwise_tukey(dv='score', between='group', data=df)
    print(posthoc)  # 包含每对比较的Hedges' g

带诊断的线性回归

import statsmodels.api as sm
from assumption_checks import check_regression_diagnostics

X = sm.add_constant(X_predictors)  # 添加截距
model = sm.OLS(y, X).fit()
print(model.summary())

# 4面板残差图 + 夏皮罗-威尔克检验、Breusch-Pagan检验、Durbin-Watson检验、VIF
diag = check_regression_diagnostics(model)
print(diag['interpretation'])
print(diag['vif'])

# 如果标记出异方差性，改为报告稳健标准误
robust = model.get_robustcov_results('HC3')

贝叶斯T检验

import pymc as pm
import arviz as az
import numpy as np

with pm.Model() as model:
    # 先验
    mu1 = pm.Normal('mu_group1', mu=0, sigma=10)
    mu2 = pm.Normal('mu_group2', mu=0, sigma=10)
    sigma = pm.HalfNormal('sigma', sigma=10)

    # 似然
    y1 = pm.Normal('y1', mu=mu1, sigma=sigma, observed=group_a)
    y2 = pm.Normal('y2', mu=mu2, sigma=sigma, observed=group_b)

    # 派生量
    diff = pm.Deterministic('difference', mu1 - mu2)

    trace = pm.sample(2000, tune=1000)

# ArviZ 1.x默认使用89%区间；报告时需显式请求95%
print(az.summary(trace, var_names=['difference'], ci_prob=0.95))

# 直接的概率陈述（这正是单侧问题所转化成的结果）
prob_greater = np.mean(trace.posterior['difference'].values > 0)
print(f"P(mu1 > mu2 | data) = {prob_greater:.3f}")

# ArviZ 1.x已移除az.plot_posterior；改用plot_dist（在0.x版本中plot_posterior仍可用）
az.plot_dist(trace, var_names=['difference'], ci_prob=0.95)

先验的尺度应与数据相匹配（例如，sigma=10适合标准差接近10的结局变量；可用观测到的标准差作为参考），并在报告中说明所用的先验。

效应量

效应量量化效应大小；p值只能说明效应是否存在。 每次检验都应报告一个效应量。完整指南参见references/effect_sizes_and_power.md。

快速参考：常见效应量

| 检验 | 效应量 | 小 | 中 | 大 |

|------|--------|----|----|----|

| T检验 | Cohen's d | 0.20 | 0.50 | 0.80 |

| 方差分析 | η²_p | 0.01 | 0.06 | 0.14 |

| 相关性 | r | 0.10 | 0.30 | 0.50 |

| 回归 | R² | 0.02 | 0.13 | 0.26 |

| 卡方检验 | Cramér's V | 0.07 | 0.21 | 0.35 |

这些基准是约定俗成的惯例，不是铁律——一个"小"效应可能极为重要（如药物副作用），而一个"大"效应也可能无关紧要。应结合具体情境解读。

计算效应量

Pingouin会在检验结果中一并返回效应量（pg.ttest的cohen_d、pg.anova的np2、pg.pairwise_tukey的hedges；pg.corr返回的r本身就是一种效应量）。

效应量的置信区间

应报告效应量的置信区间以体现其精度。使用pg.compute_esci（注意：pg.compute_effsize_from_t只返回点估计——它不会返回置信区间）：

import pingouin as pg

d = pg.compute_effsize(group_a, group_b, eftype='cohen')
ci_lower, ci_upper = pg.compute_esci(stat=d, nx=len(group_a), ny=len(group_b),
                                     eftype='cohen', confidence=0.95)
print(f"d = {d:.2f}, 95% CI [{ci_lower:.2f}, {ci_upper:.2f}]")

功效分析

先验功效分析（研究规划）

在数据收集前确定所需的样本量：

from statsmodels.stats.power import tt_ind_solve_power, FTestAnovaPower

# T检验：检测d = 0.5需要每组多少n？
n_required = tt_ind_solve_power(
    effect_size=0.5,
    alpha=0.05,
    power=0.80,
    ratio=1.0,
    alternative='two-sided'
)
print(f"每组所需n: {n_required:.0f}")

# 单因素方差分析：检测Cohen's f = 0.25需要多少n？
# 注意：参数名是k_groups；effect_size是Cohen's f（f = sqrt(eta2/(1-eta2)));
# 且solve_power返回的是总样本量，而不是每组的n。
import math
anova_power = FTestAnovaPower()
n_total = anova_power.solve_power(
    effect_size=0.25,
    k_groups=3,
    alpha=0.05,
    power=0.80
)
print(f"所需总样本量: {math.ceil(n_total)}（每组约{math.ceil(n_total / 3)}）")

敏感性分析（研究后）

确定该研究能够检测到的效应量：

# 每组n=50时，在80%功效下我们可以检测到什么效应？
detectable_d = tt_ind_solve_power(
    effect_size=None,  # 求解这个
    nobs1=50,
    alpha=0.05,
    power=0.80,
    ratio=1.0,
    alternative='two-sided'
)
print(f"研究可以检测到d >= {detectable_d:.2f}")

注意：事后的"观测功效"（根据观测到的效应量计算功效）是循环论证且具有误导性的——它其实是p值的确定性函数。如果研究已经完成，有人询问功效问题，应改为进行敏感性分析。

有关详细指导，请参见references/effect_sizes_and_power.md。

报告结果

遵循references/reporting_standards.md中的APA风格指南。每份报告都需要包含：

1. 描述性统计：所有组/变量的M、SD、n
2. 检验统计量：检验名称、统计量、df、精确p值（写作p = .034，而不是p < .05；只有p值低于.001时才使用p < .001）
3. 效应量：附带置信区间
4. 假设检查：进行了哪些检验、结果如何、采取了什么行动
5. 所有计划的分析：包括不显著的结果——省略它们就是选择性报告

示例报告模板

独立T检验

组A（n = 48, M = 75.2, SD = 8.5）的得分显著高于
组B（n = 52, M = 68.3, SD = 9.2），t(98) = 3.82, p < .001, d = 0.77,
95% CI [0.36, 1.18]，双侧检验。正态性假设（夏皮罗-威尔克检验：
组A W = 0.97, p = .18；组B W = 0.96, p = .12）和方差齐性假设
（列文检验F(1, 98) = 1.23, p = .27）均得到满足。

单因素方差分析

单因素方差分析显示治疗条件对测试得分有显著主效应，
F(2, 147) = 8.45, p < .001, η²_p = .10。使用Tukey's HSD的事后
比较表明，条件A（M = 78.2, SD = 7.3）的得分显著高于
条件B（M = 71.5, SD = 8.1, p = .002, d = 0.87）和
条件C（M = 70.1, SD = 7.9, p < .001, d = 1.07）。
条件B和C之间无显著差异（p = .52, d = 0.18）。

多元回归

进行多元线性回归以从学习时间、先前GPA和出勤率预测考试分数。
整体模型显著，F(3, 146) = 45.2, p < .001, R² = .48, 调整R² = .47。
学习时间（B = 1.80, SE = 0.31, β = .35, t = 5.78, p < .001, 95% CI [1.18, 2.42]）
和先前GPA（B = 8.52, SE = 1.95, β = .28, t = 4.37, p < .001, 95% CI [4.66, 12.38]）
是显著预测因子，而出勤率不是（B = 0.15, SE = 0.12, β = .08, t = 1.25, p = .21, 95% CI [-0.09, 0.39]）。
多重共线性不是问题（所有VIF < 1.5）。

贝叶斯分析

使用弱信息先验（组均值为Normal(0, 10)）进行贝叶斯独立样本t检验。
后验分布表明组A得分高于组B
（M_diff = 6.8, 95% 可信区间 [3.2, 10.4]），
组A均值超过组B均值的后验概率为99.8%。
收敛诊断令人满意（所有R-hat < 1.01, ESS > 1000）。

如果使用了非参数检验，应报告中位数而不是均值、报告U/W/H统计量，以及基于秩的效应量（例如秩双列相关系数，pg.mwu会以RBC列返回该值）。

贝叶斯统计

在以下情况下考虑使用贝叶斯方法：

• 有先验信息需要纳入分析
• 想要关于假设的直接概率陈述（"该效应有95%的概率落在此区间内"）
• 样本量较小，或数据收集是序贯进行的（无需为可选停止校正）
• 需要量化*支持*零假设的证据
• 模型较为复杂（层次结构、缺失数据）

有关先验设定、贝叶斯因子、可信区间、层次模型以及收敛检查（R-hat < 1.01、充足的ESS、后验预测检查），请参见references/bayesian_statistics.md。

捆绑资源

参考文档（references/）

• test_selection_guide.md：涵盖组间比较、关系分析、计数数据、时间-事件分析、一致性/信度和分类分析的决策树
• assumptions_and_diagnostics.md：关于检查和处理假设违反情况的详细指导
• effect_sizes_and_power.md：计算、解释和报告效应量；功效分析
• bayesian_statistics.md：先验、贝叶斯因子、可信区间、层次模型、诊断
• reporting_standards.md：带实例的APA风格报告指南

脚本（scripts/）

• assumption_checks.py：带可视化的自动假设检查
• comprehensive_assumption_check()：一次调用完成异常值 + 正态性 + 方差齐性检验
• check_normality()、check_normality_per_group()：带Q-Q图的夏皮罗-威尔克检验
• check_homogeneity_of_variance()：带箱线图的列文检验
• check_regression_diagnostics()：针对已拟合OLS模型的4面板残差图 + 夏皮罗-威尔克检验、Breusch-Pagan检验、Durbin-Watson检验、VIF
• check_linearity()、detect_outliers()

统计诚信

以下是让一份分析经得起推敲的实践。它们之所以重要，是因为最常见的统计错误并非计算失误——而是悄悄的"灵活性"（不断尝试直到得出想要的结果）和选择性报告。

1. 区分验证性分析和探索性分析。 在运行分析之前就说明计划采用的分析方法；将过程中发现的任何结果标注为探索性发现。
2. 不要为了追求显著性而反复尝试。 如果计划中的检验结果不显著，那就是结果本身。反复尝试不同的检验方法、子组或异常值剔除方案直到p < .05，会使该p值失去意义。
3. 在运行一系列检验时，应校正多重比较（事后方差分析用Tukey HSD；其他系列检验用Holm或Benjamini-Hochberg FDR校正），并说明使用了哪种校正方法。
4. 不显著的结果不等于没有效应的证据。 在样本量较小的情况下，研究可能只是功效不足——应运行敏感性分析，或使用贝叶斯分析/等效性检验来真正量化对零假设的支持程度。
5. 统计显著性不等于实际重要性。 在样本量较大时，微不足道的效应也能达到p < .001。解读时应以效应量为主导。
6. 在删除行之前先理解缺失数据的机制。 列表式删除（listwise deletion）只有在数据完全随机缺失（MCAR）时才是安全的；否则应考虑多重插补，并说明所采取的处理方式。
7. 保证可复现性。 设置随机种子，对基于模拟的方法报告所用的库版本，并将分析保存为可运行的脚本。