1. 相关性分析真正在回答的问题相关性分析这个词被用得太泛了。做数据的人嘴里说的相关有时候是统计学意义上的相关系数有时候只是这两列看起来有关系还有时候是因果关系的委婉说法。我见过不少分析报告拿一个 Pearson 系数就下结论说A 越高 B 越大但实际上 A 是个性别字段被编码成 0/1B 是用户消费金额——这种相关系数算出来的数字在数学上没毛病在业务解释上基本是一场灾难。问题就出在变量类型上。相关性分析的本质是衡量两个变量共同变化的程度和方向而共同变化这四个字对不同类型的变量有不同的数学定义。连续变量之间的共同变化可以用协方差刻画分类变量之间的共同变化只能靠频数分布去比对而分类变量和连续变量之间的共同变化既没有协方差可言、也没有频数表可列必须借助组间差异这个中间桥梁来翻译。这就是为什么分类变量与连续变量的相关性分析是最容易翻车的一块它不像连续对连续那样有唯一的标准答案二分类、多分类、有序分类各自对应不同的方法选错了数字还能跑出来只是结论全错。这篇东西就是想把这块讲透。适合谁看如果你做过用户画像、AB 实验分析、问卷数据清洗、风控建模前的特征筛查或者你只是单纯想知道性别和收入到底有没有关系、关系多大该怎么严谨地回答那下面的内容应该能直接用上。我会按连续对连续、分类对分类、分类对连续三条线走重点压在分类对连续上因为那是方法最多、坑也最多的地方最后给一套能跑的代码和一张排查表。1.1 变量类型决定了你能问什么问题先说一个容易被忽略的前提变量类型不是数据自带的属性而是分析者根据业务含义赋予的。同样一列用户等级取值 1 到 5你可以把它当有序分类等级高低有意义但 1 到 2 的距离不等于 2 到 3 的距离也可以当连续变量假设间隔均匀甚至当无序分类只关心属于哪一档。这三种身份会导向三种完全不同的分析方法结论也可能不一样。这个认识很重要。我见过太多人拿到数据第一件事就是df.corr()然后对着一张全是数字的热力图发呆完全没意识到里面有些列是标签编码后的省市区、有些是独热编码的残渣、有些是 0/1 的布尔值。corr()默认跑 Pearson它对输入的唯一要求是能转成数值转完之后它不管语义一律按等距变量处理。分类变量被编码成整数之后Pearson 算出来的东西没有统计学意义——它衡量的是编码值的线性关系而不是类别与数值的关系。省市区编码里北京1、上海2算出来的相关系数你跟谁都没法解释。所以相关性分析的第一步永远是分类把每一列贴上连续有序分类无序分类三个标签之一然后按标签组合去查方法表。方法表长这样变量 A 类型变量 B 类型推荐方法效应量指标对称性连续连续Pearson / Spearmanr / ρ对称无序分类无序分类卡方检验 / Fisher 精确Cramérs V、φ对称二分类连续点二列相关 / 独立样本 t 检验r_pb、Cohens d对称多分类连续单因素方差分析 / Kruskal-Wallisη²、ε²不对称有序分类连续Spearman / Jonckheere-Terpstraρ、趋势统计量对称有序分类有序分类Spearman / Kendalls τ-bρ、τ对称表里对称性那一列值得单独说一句。Pearson、Spearman、Cramérs V 都是对称的A 和 B 换个位置数字不变。但方差分析导出的 η² 不是——它算的是分类变量能解释连续变量多少变异反过来问连续变量能解释分类变量多少就没法用 η² 表达了。这种不对称在写报告时要注意措辞别写成两者互相影响那已经滑向因果论了。1.2 分类与连续之间的相关为什么总被算成假的我总结过三类高频的假相关几乎每个项目都会撞上一个。第一类是编码当连续。把无序分类做 LabelEncoder 之后直接进corr()。比如城市编码、渠道 ID、商品类目 ID这些都只是编号编号大小没有任何顺序含义。算出来的 Pearson 系数纯属数字巧合样本量一大还能显著特别能唬人。判断方法很简单看这一列的取值是否等距、是否有序、0 值是否代表没有。任何一个不成立就别用 Pearson。第二类是把连续变量硬切成组再算相关。这属于信息损失型错误。有人为了方便看把年龄切成青年/中年/老年三组然后算这三组和收入的关系。问题在于切点是你自己定的换个切点结论就可能变而且分组会大幅损失变异信息导致检验功效下降。原本连续对连续也许能看出 0.2 的相关切完组之后组间差异不显著你就误判成没关系。这个错误的方向是保守的但保守也意味着漏掉真信号。第三类是忽略分布形态。多分类对连续标准做法是方差分析而方差分析对正态性和方差齐性有要求同时对极端值极其敏感。一个组里混进两个异常大的值均值被拉高F 值就上去了你得到组间差异显著回头查数据发现是两个录错的 999999。所以做这类分析之前箱线图必须画Levene 检验必须跑扭不过去的就换 Welch ANOVA 或者 Kruskal-Wallis。这三类错误有个共同点都能顺利跑出数字都不报错都看起来很像那么回事。这就是为什么我在团队里立了个规矩凡是涉及分类变量的相关性分析结果必须附上方法选择说明——用了什么方法、为什么用它、前提检查过没有。没有这段说明的结果一律打回。2. 变量类型判定与数据前置体检方法选对了一半剩下一半靠数据本身干不干净。相关性分析对数据质量的要求比很多人想象的高缺失值处理方式不同相关系数能差出一大截异常值不处理Pearson 能被单个点带跑类别频数太悬殊卡方和方差分析都会失真。这一节讲的就是上桌前必须做的事。2.1 连续、有序、名义三种身份的边界判定的核心问题是三个取值能不能比大小、差值是否等距、零点是否有实在含义。三个都是是就是连续变量收入、时长、温度、金额。能比大小但差值不等距是有序分类满意度五级、学历层次、信用等级。连大小都不能比就是无序分类城市、渠道、设备型号、是否点击。边界上有几个常被搞混的例子值得单独拎出来李克特量表。问卷里常见的非常不满意到非常满意五分制绝大多数场合应该当有序分类处理。虽然实践中很多人直接当连续变量算均值但那需要额外假设相邻选项之间的心理距离相等这个假设在题目设计良好时可以勉强成立设计粗糙时就是瞎猜。稳妥起见做相关分析时用 Spearman 而不是 Pearson两种做法的结论一致时再放心用均值。金额分段。业务报表里常见的0-100 元、100-500 元、500-1000 元、1000 元以上这种分箱本质上是把连续变量切成有序分类。切完之后不该再当连续用因为箱宽不一致等距假设直接崩了。布尔值。0/1 变量既可以看作二分类也可以看作特殊的连续变量。这个双重身份很有用——正是因为二分类可以安全地编码成 0/1 而不损失信息点二列相关才能在数学上等价于 Pearson 相关系数。但前提是它真的只有两个取值而不是金额为 0 或大于 0这种被人工二分化的产物后者属于人为切分得用二列相关后面会讲。ID 类字段。订单号、用户 ID 这种看着是数值实际是名义分类而且类别数等于样本量任何基于频数的方法都会崩。这类字段在相关性分析里应该直接排除不要试图分析。判定完之后建议在数据字典里固化下来写成枚举而不是靠记忆。项目一多人一多同一个字段被两个人当两种类型用的情况太常见了最后两份报告的结论对不上排查成本极高。2.2 上桌前先做的四件事第一件缺失值体检。相关性分析对缺失的处理有三种常见策略——按列删除、按对删除、插补。pandas.corr()默认是成对删除pairwise这意味着不同格子的相关系数是基于不同样本量算出来的严格来说不可比。样本缺失有结构性的时候更麻烦如果高收入人群普遍不填收入那收入与消费的相关实际上只在中低收入群体里算出来的解释时得打个大问号。我的习惯是先看缺失比例和缺失是否与其他变量相关缺失超过 20% 的列直接标记为谨慎使用。第二件异常值体检。连续变量用箱线图或者 1.5 倍 IQR 规则扫一遍重点看不合理的大值比如年龄 200、时长负数。异常值对 Pearson 的影响是压倒性的对 Spearman 影响小很多对基于秩的方法几乎无影响。所以如果你发现 Pearson 和 Spearman 结论不一致八成是异常值在作祟这时候要回到原始数据去看而不是在两个系数之间挑一个好看的去汇报。第三件类别频数体检。分类变量看每一类的样本量。卡方检验的经验规则是期望频数不低于 5 的格子要占 80% 以上方差分析虽然对组样本量的不平衡有一定容忍度但一组只有 3 个样本、另一组有 300 个均值的稳定性完全不在一个量级F 检验的结果就不可靠了。遇到小类别要么合并注意合并要有业务依据别把高风险和中风险合并了要么改用精确检验或者非参数方法。第四件分布形态体检。连续变量看偏度和峰度或者直接上直方图。严重右偏的收入、时长类变量取对数之后往往更接近正态这时候对数变换后的 Pearson 和原始的 Spearman 结论通常一致可以互相印证。方差分析前要跑正态性检验Shapiro 或 Anderson-Darling和方差齐性检验Levene样本量大时正态性检验会过度敏感这时候看 Q-Q 图比看 p 值靠谱。这四件事做完其实你已经对数据有了基本判断该用哪个方法、结果能信几分心里大概有数了。跳过这步直接算后面返工的概率接近百分之百。3. 连续对连续拿来做对照的基线分类对连续的方法很多是从连续对连续的方法里变形出来的点二列相关是 Pearson 的特例方差分析本质上是组间方差与组内方差的比值η² 是 r² 的多组推广。所以先把连续对连续这条基线捋清楚后面的推导会顺畅很多。已经熟的人可以跳过这一节但下面关于置信区间和样本量的部分我建议还是扫一眼因为那是分类对连续里同样存在的陷阱。3.1 Pearson 与 Spearman 的取舍逻辑Pearson 相关系数衡量的是线性关系公式是协方差除以两个标准差的乘积取值在 -1 到 1 之间。它对线性这两个字非常敏感如果两个变量是完美的 U 型关系Pearson 会算出接近 0 的结果告诉你没关系这显然是错的。所以 Pearson 只能在肉眼确认散点图大致呈线性分布之后用。Spearman 是 Pearson 在秩上的版本把两个变量各自排序用秩次代替原值算 Pearson。它的好处是只依赖顺序不依赖具体数值所以对单调非线性关系比如对数关系也能捕捉到同时对异常值免疫。代价是它丢掉了距离信息两个变量都是等距变量、关系确实线性时Spearman 的功效略低于 Pearson。实际项目里的取舍原则我一般这么定散点图线性、无明显异常值、双变量近似正态 → Pearson单调但非线性、有异常值、变量是有序分类 → Spearman样本量小小于 30→ Spearman 更稳因为 Pearson 对正态假设的依赖在这个样本量下很明显两边都要算结论一致就放心报 Pearson不一致就去看散点图找原因两个都算这个做法看起来偷懒其实是最实用的诊断手段。两者的差值本身就是信息差了 0.1 以上说明数据里要么有异常值、要么关系非线性值得进一步查。顺便说一句 Kendalls τ。它和 Spearman 都是基于秩的但计算逻辑不同Spearman 是秩的 PearsonKendall 是一致对减不一致对的比例。Kendall 在小样本和大量并列值ties的情况下表现更稳解释也更直观——τ 0.3 大致意味着随机抽两对观测同向的概率比反向高 30%。缺点是计算复杂度高大样本时慢。样本量上万的数据用 Spearman 就够了。3.2 显著性不等于可靠区间估计与样本量这部分是我最想强调的。p 0.05 只说明在假设没有关系的前提下观测到这个系数的概率很小它不说明关系有多强。样本量一大r 0.05 也能显著样本量一小r 0.6 也可能不显著。拿 p 值当有没有关系的判据是相关性分析里最常见的误用。正确的做法是同时看三样东西系数大小、置信区间、样本量。相关系数的置信区间可以用 Fisher z 变换来算。先把 r 转成 z arctanh(r)z 的标准误是 1/√(n-3)据此算区间再反变换回 r。代码很简单import numpy as np from scipy import stats def corr_ci(r, n, alpha0.05): z np.arctanh(r) se 1 / np.sqrt(n - 3) crit stats.norm.ppf(1 - alpha / 2) lo, hi np.tanh([z - crit * se, z crit * se]) return lo, hi print(corr_ci(0.15, 2000)) # 窄区间说明估计精细 print(corr_ci(0.15, 40)) # 宽区间跨过 0说明样本不足同样 r 0.152000 个样本时区间大概是 [0.11, 0.19]40 个样本时区间可能从 -0.17 到 0.44直接跨零。这两种情况下r 0.15这句话的可信度天差地别但很多报告只写一个数字。效应量的经验参照不同领域差异很大仅作粗略参考|r| 在 0.1 附近算小0.3 中等0.5 以上算大。社会科学里 0.2 已经值得讨论物理实验里 0.99 才算合格。所以别照搬别人的阈值要结合自己领域的常识。还有一个常被忽略的点相关系数会被测量误差稀释。如果两个变量都有随机测量误差观测到的相关系数会系统性地低于真实相关衰减偏误。问卷数据尤其明显一道题的测量信度只有 0.7两个变量各有一半误差观测相关就会被压低不少。所以拿到一个偏低的系数别急着下没关系的结论先看看变量本身测得准不准。4. 分类对分类卡方与它的效应量分类对分类的情况在业务里非常常见性别与是否流失、渠道与是否转化、会员等级与复购档位。这块的方法相对统一但同样有一个高频错误——只报 p 值不报效应量。4.1 卡方检验成立的前提卡方检验的逻辑是假设两个分类变量独立按行合计和列合计推算每个格子的期望频数然后比较观测频数与期望频数的偏离程度。偏离越大越不像独立。统计量是 Σ(观测-期望)²/期望自由度是 (行数-1)×(列数-1)。用之前必须检查的三条期望频数。这是硬条件。任一格子的期望频数小于 1或者超过 20% 的格子期望频数小于 5卡方分布的近似就不可靠了p 值会失真。处理方式合并类别合并要有业务逻辑不能随便合或者改用 Fisher 精确检验2×2 表完全适用更大的表可以用 Fisher-Freeman-Halton计算量大但现代工具都能跑。样本独立性。每个观测只能落在一个格子里不能重复计数。同一个用户既算新客又算老客的话卡方就废了。这个在宽表拼接的时候特别容易出错多表 join 之后行数膨胀重复用户被数了两次。类别定义一致。两个变量的分类标准要在时间上和口径上对齐。这个不是统计问题而是数据治理问题但影响更大——去年按订单维度统计的复购今年按用户维度统计两张表放一起做卡方结果没有意义。4.2 从 p 值走到 Cramérs V卡方统计量本身不能当效应量用因为它的数值随样本量线性增长同样的分布偏离样本量翻倍卡方值翻倍。所以需要一个归一化的指标这就是 Cramérs VV sqrt(χ² / (n × (min(r, c) - 1)))取值 0 到 10 表示完全独立。2×2 表时 V 就等于 φ 系数这时它的符号可以表示方向正相关/负相关但 Cramérs V 在更大的表里没有方向概念只有强度。这里有个细节要注意Cramérs V 在小样本、大表的情况下会有向上的偏误即真实无关的两个变量也可能算出 0.1 到 0.2 的 V。有个 Bergsma 校正公式可以修正实践中如果表大于 3×3 且样本量不足 200我建议用校正版def cramers_v_corrected(x, y): table pd.crosstab(x, y) chi2 stats.chi2_contingency(table, correctionFalse)[0] n table.values.sum() r, k table.shape phi2 chi2 / n phi2corr max(0, phi2 - (k - 1) * (r - 1) / (n - 1)) rcorr r - (r - 1) ** 2 / (n - 1) kcorr k - (k - 1) ** 2 / (n - 1) denom min(kcorr - 1, rcorr - 1) return np.sqrt(phi2corr / denom) if denom 0 else np.nan效应量大小的经验阈值2×2 表时 φ 的 0.1/0.3/0.5 对应小/中/大自由度更大时要按自由度调整粗略的做法是把 V 除以 √(min(r,c)-1) 再比对。这块的争议一直存在不同教材给的阈值不一样我倾向的做法是不要死守阈值而是拿同类研究或历史数据做参照——本次 V 0.12去年同口径是 0.11基本稳定这种表述比V 0.12 属于小效应有用得多。还有一点必须强调卡方显著不等于有关系更不等于有关系强度。n 10000 的时候V 0.03 都能显著。这种统计显著但业务无感的结论在流量数据里遍地都是。所以每次报卡方我都会强制附上 V 和交叉表的百分比分布让看图的人自己判断差异有没有业务意义。5. 分类对连续主线方法拆解进入正题。分类变量与连续变量的相关性分析核心思路是把类别翻译成组把相关翻译成组间差异相对于组内波动的比例。翻译得对不对取决于分类变量有几个水平、这些水平有没有顺序。下面按二分类、多分类、有序分类三条支线拆最后单独讲人为二分化的坑。5.1 二分类对连续点二列、t 检验与 Cohens d二分类变量性别、是否会员、是否曝光和连续变量消费金额、停留时长的关系最标准的指标是点二列相关系数r_pb。它的定义和 Pearson 完全等价——只要把二分类编码成 0/1两者数值一模一样。之所以单独立一个名字是因为它的平方有明确的解释r_pb² 等于组间平方和占总平方和的比例也就是这一组 0/1 标签能解释连续变量多少变异。计算方式有两种结果相同import numpy as np from scipy import stats # 方式一直接调 scipy r_pb, p_pb stats.pointbiserialr(df[is_member], df[amount]) # 方式二手算顺便验证公式 n1 (df[is_member] 1).sum() n0 (df[is_member] 0).sum() m1 df.loc[df[is_member] 1, amount].mean() m0 df.loc[df[is_member] 0, amount].mean() s df[amount].std(ddof1) n n1 n0 r_manual (m1 - m0) / s * np.sqrt(n0 * n1 / (n * (n - 1)))公式里的 (m1 - m0)/s 就是标准化均值差后面的根号项是一个和组样本量有关的系数。这告诉我们一个重要事实r_pb 的大小同时受均值差有多大和两组样本量是否均衡影响。同样差 0.5 个标准差两组各 50 人的 r_pb 会明显大于 90 人 vs 10 人的情况。所以两组人数悬殊时r_pb 会偏小别误读成关系弱。与之对应的、不受样本量比例影响的指标是Cohens dd (m1 - m0) / s_pooledd 只反映差了几个标准差0.2/0.5/0.8 对应小/中/大。在 AB 实验报告里d 比 r_pb 更常用因为它稳定、可比、和样本量解耦。我的习惯是两个都报d 给业务看差距大小r_pb² 给建模看解释力。t 检验和点二列相关的关系也得说清楚。独立样本 t 检验检验的是两组均值是否相等它和 r_pb 的显著性检验在数学上等价t² (n-2) r_pb² / (1 - r_pb²)。换句话说同一个问题t 检验和点二列相关给出完全一致的 p 值只是一个从均值差角度讲一个从共变角度讲。别把两者当成两种独立证据去交叉验证它们是同一个东西的两副面孔。实操上要注意 t 检验的方差齐性。默认的 Student t 检验假设两组方差相等不满足时会偏。直接用 Welch t 检验equal_varFalse它在方差不齐时表现更好方差齐时损失的功效很小所以我现在基本默认用 Welch。样本量很小或者分布严重偏斜时改用 Mann-Whitney U 检验配套的效应量用秩二列相关u, p stats.mannwhitneyu(g1, g0, alternativetwo-sided) r_rb 1 - 2 * u / (len(g1) * len(g0)) # 取值 -1 到 15.2 多分类对连续方差分析、eta 方与 Welch 兜底三个及以上类别会员等级、城市线级、渠道来源对连续变量标准方案是单因素方差分析。它的核心思想是分解方差总变异 组间变异 组内变异。如果组间变异相对于组内变异足够大就认为各组均值不全相等。F 统计量 (组间均方) / (组内均方)自由度是 (k-1, n-k)。这里 k 是组数n 是总样本量。F 显著只说明至少有一组和其他组不同不说明是哪一组也不说明差异有多大——所以必须做事后多重比较而且必须校正。效应量用η²eta squaredη² SS_between / SS_total解释是分类变量解释了连续变量百分之多少的变异取值 0 到 1。它和 r² 是同一个东西在多组情形下的推广二分类时 η² 就等于 r_pb²。这个是很多人不知道的细节但很有用它意味着你可以把不同方法算出来的效应量放在同一把尺子上比较。η² 有个已知问题它是有偏的样本量小时会高估。所以还有个ω²omega squared和偏 η²partial eta squared。偏 η² 是统计软件默认输出的那个比如某些工具的 ANOVA 表但它在单因素设计里和 η² 数值一样多因素设计里才不同。报告时写清楚用的是哪个别混。手算 η² 很简单比调包还快grand_mean df[amount].mean() ss_between sum(len(g) * (g[amount].mean() - grand_mean) ** 2 for _, g in df.groupby(level)) ss_total ((df[amount] - grand_mean) ** 2).sum() eta2 ss_between / ss_total f_stat, p_val stats.f_oneway(*[g[amount] for _, g in df.groupby(level)])前提检查两条正态性每组残差近似正态和方差齐性Levene 检验 p 0.05。方差不齐时用Welch ANOVA它不假设方差相等自由度用 Welch-Satterthwaite 近似import pingouin as pg pg.welch_anova(datadf, dvamount, betweenlevel)严重偏态、异常值多、或者有序分类当无序处理时换Kruskal-Wallis检验。它是方差分析的非参数版本比较的是各组的秩均值。显著之后的事后比较用 Dunn 检验配 Bonferroni 或 BH 校正效应量可以用ε²epsilon squaredε² (H - k 1) / (n - k)其中 H 是 Kruskal-Wallis 统计量k 是组数n 是总样本量。ε² 的解释和 η² 类似但通常略小。事后比较这一步千万别省也千万别用 t 检验连打好几次还不校正。k 组两两比较有 k(k-1)/2 次5 组就是 10 次不做校正的话假阳性率会从 5% 飙到 40% 左右。校正方法用 Bonferroni保守比较次数少时够用或者 Benjamini-Hochberg控制 FDR比较次数多时更合适。5.3 有序分类对连续别把等级当哑变量有序分类是三种类型里最尴尬的它比无序分类多了一层顺序信息比连续变量少了一层距离信息。处理方式取决于你想问什么问题。如果只想问有没有差异把它当无序分类走方差分析或 Kruskal-Wallis。这样做不损失正确性但浪费了顺序信息功效会低一些。如果想问有没有单调趋势那就要用能利用顺序的方法。最直接的是把等级编码成 1、2、3、4、5然后和连续变量算Spearman 秩相关。这个做法在等级数不少于 5 且分布大致均匀时表现不错等级只有 3 档时会比较粗糙但方向性判断仍然可用。更严谨的做法是Jonckheere-Terpstra 趋势检验它专门检验随着等级升高连续变量的分布是否有序移动不要求等距假设功效比 Kruskal-Wallis 高。这个检验在常规工具里不一定有R 里有clinfun::jonckheere.testPython 里可以用scikit-posthocs或者自己按公式实现。绝对要避免的做法是把有序分类当无序分类拆成哑变量然后逐个看系数。这样得到的是一堆等级 2 比等级 1 高多少等级 3 比等级 1 高多少的比较完全丢掉了等级越高越好这个整体趋势解释起来支离破碎。如果确实要进回归模型用正交多项式对比或者直接当连续变量处理前提是等距假设能接受都比拆哑变量强。有一次做用户等级和客单价的分析一开始按五个等级拆哑变量跑 ANOVAF 检验显著但事后比较的结果很乱看不出规律。后来换成 Spearmanρ 0.31p 值极小一句话就说清了等级越高客单价越高。同一个数据换个方法可解释性天差地别。5.4 手工二分化连续变量的代价二列相关这个坑值得单独开一节因为它太容易被当成正常操作。场景是这样的你有一个连续变量比如用户活跃度得分业务上想按某个阈值切成高活跃/低活跃然后分析它和另一个连续变量消费金额的关系。切完之后你算了个点二列相关得到 r 0.28。这个数字是低估的。原因是二分化丢掉了组内的全部变异信息。原本活跃度从 10 分到 90 分是一个连续谱切完之后只剩高/低两个标签中间那么多信息全没了。统计上这叫人为二分化的信息损失会导致相关系数系统性偏小除非切点在均值处且分布正态。如果你手上就是二分后的数据对应的正确指标是二列相关系数biserial correlation它估计的是假如那个连续变量没被切原本的相关会是多少。公式是r_bis r_pb × sqrt(p × q) / φ(t)其中 p 和 q 是两组的比例φ(t) 是标准正态分布在切分点处的密度值切分点 t 是使正态分布上侧面积为 p 的那个 z 值。算出来 r_bis 通常会比 r_pb 大一些因为它在补偿信息损失。具体怎么算from scipy.stats import norm def biserial(r_pb, p): q 1 - p t norm.ppf(q) # 切分点 density norm.pdf(t) # 该点处的正态密度 return r_pb * np.sqrt(p * q) / density p_high (df[active_flag] 1).mean() print(biserial(0.28, p_high))注意二列相关的假设很强它要求被二分的那个变量本身服从正态分布且切分是人为的不是天然的二分变量。如果活跃度本身就是双峰的或者高/低是天然存在的两类用户那就该用点二列不该用二列硬套会得到错误结果。我的实际建议是能不切就不切。业务上需要高/低这种标签可以在报告层做分组展示但分析层始终保留连续形式。真要切切点要有依据比如业务阈值、分位数而且要意识到所有基于二分数据的相关系数都是下界估计。6. 完整实操一套能跑的相关性矩阵前面讲的是分散的方法实际项目里真正有用的是输入一张表自动按变量类型分发到正确方法输出统一的效应量和 p 值。这一节给一套能直接抄的代码。6.1 环境准备与模拟数据依赖三个库pandas 做数据处理scipy 做核心检验pingouin 补充 Welch ANOVA 和偏 η² 这类 scipy 里没有的实现。pip install pandas numpy scipy pingouin scikit-posthocs造一份覆盖三种变量类型的模拟数据。为了让后面的结果可验证我按已知的组间差异生成数据import numpy as np import pandas as pd from scipy import stats import pingouin as pg rng np.random.default_rng(20240607) n 480 # 无序分类三个渠道 channel rng.choice([自然, 投放, 推荐], sizen, p[0.5, 0.3, 0.2]) # 有序分类五个会员等级 level rng.choice([1, 2, 3, 4, 5], sizen, p[0.3, 0.25, 0.2, 0.15, 0.1]) # 二分类是否会员 is_member (level 3).astype(int) # 连续变量消费金额受渠道和等级共同影响加噪声 base np.where(channel 自然, 120, np.where(channel 投放, 145, 165)) amount base level * 12 rng.normal(0, 45, n) amount np.clip(amount, 0, None) df pd.DataFrame({ channel: channel, level: level, is_member: is_member, amount: amount, visits: rng.normal(8, 3, n) level * 0.8 })这份数据里我埋了两条已知的关系渠道影响金额、等级影响金额后面就看方法能不能把它们捞出来。6.2 逐类方法的代码与结果解读先跑二分类对连续# 点二列相关 r_pb, p_pb stats.pointbiserialr(df[is_member], df[amount]) print(f点二列 r {r_pb:.3f}, p {p_pb:.4g}) # Cohens dWelch g1 df.loc[df[is_member] 1, amount] g0 df.loc[df[is_member] 0, amount] sp np.sqrt(((len(g1)-1)*g1.var(ddof1) (len(g0)-1)*g0.var(ddof1)) / (len(g1) len(g0) - 2)) d (g1.mean() - g0.mean()) / sp print(fCohens d {d:.3f}) # Welch t 检验 t, p_t stats.ttest_ind(g1, g0, equal_varFalse) print(fWelch t {t:.3f}, p {p_t:.4g})跑出来应该是 d 在 0.9 到 1.2 之间r_pb 在 0.4 上下。这里能验证前面说的等价性t² 和 r_pb 换算出来的值应该严丝合缝对得上verificationt² (n-2)r²/(1-r²)。再跑多分类对连续# 方差分析 groups [g[amount].values for _, g in df.groupby(channel)] f_stat, p_anova stats.f_oneway(*groups) # eta 方 grand df[amount].mean() ss_between sum(len(g) * (g[amount].mean() - grand) ** 2 for _, g in df.groupby(channel)) ss_total ((df[amount] - grand) ** 2).sum() eta2 ss_between / ss_total print(fF {f_stat:.3f}, p {p_anova:.4g}, eta2 {eta2:.3f}) # 方差齐性 lev stats.levene(*groups) print(fLevene p {lev.pvalue:.4g}) # 事后比较Dunn BH 校正非参数路线 import scikit_posthocs as sp_ print(sp_.posthoc_dunn(df, val_colamount, group_colchannel, p_adjustfdr_bh))注意看 Levene 的结果。方差齐的话 ANOVA 没问题不齐就换pg.welch_anova。这份模拟数据里三组的标准差都设成 45理论上应该齐但抽样波动可能让 p 值在 0.05 附近晃这时候别纠结直接上 Welch 更省心。有序分类对连续rho, p_rho stats.spearmanr(df[level], df[amount]) print(fSpearman rho {rho:.3f}, p {p_rho:.4g}) # 对照当连续变量处理算 Pearson r_pear, p_pear stats.pearsonr(df[level], df[amount]) print(fPearson r {r_pear:.3f})两个数应该很接近因为我在生成数据时让等级对金额的影响是线性的这验证了等距假设成立时两种方法一致。如果真实数据里等级对金额是加速上升的比如等级 5 突然贵很多Spearman 会明显高于 Pearson那时候更要信 Spearman。6.3 组装成自动分发的相关矩阵把上面的逻辑封装成函数按变量类型自动选方法。核心是一个类型判定函数和一个分发函数def detect_type(s, max_cat10): if s.dtype object or s.dtype.name category: n_unique s.nunique() return ordinal if is_ordered(s) else nominal n_unique s.nunique() if n_unique 2: return binary if n_unique max_cat and (s s.astype(int)).all(): return ordinal # 少数几个整数档默认当有序 return continuous def pair_assoc(df, col_a, col_b): ta, tb detect_type(df[col_a]), detect_type(df[col_b]) a, b df[col_a], df[col_b] # 统一成分类在前、连续在后的顺序处理 if ta continuous and tb ! continuous: a, b, ta, tb b, a, tb, ta if ta continuous and tb continuous: r, p stats.pearsonr(a, b) return {method: pearson, effect: r, p: p} if ta in (nominal, ordinal) and tb nominal and ta nominal: table pd.crosstab(a, b) chi2, p, _, _ stats.chi2_contingency(table) v np.sqrt(chi2 / (len(a) * (min(table.shape) - 1))) return {method: cramers_v, effect: v, p: p} if tb continuous: if ta binary: r, p stats.pointbiserialr(a, b) return {method: point_biserial, effect: r, p: p} if ta ordinal: r, p stats.spearmanr(a, b) return {method: spearman, effect: r, p: p} f, p stats.f_oneway(*[g[b.name].values for _, g in df.groupby(a.name)]) grand b.mean() eta2 sum(len(g) * (g[b.name].mean() - grand) ** 2 for _, g in df.groupby(a.name)) / ((b - grand) ** 2).sum() return {method: anova_eta2, effect: eta2, p: p} return {method: unsupported, effect: np.nan, p: np.nan}这个函数不完美is_ordered需要你根据业务字典来实现类型判定也有模糊地带但它把方法论固化进了代码比每次人工判断靠谱得多。真要用的话补三件事多重比较校正循环里累积 p 值之后统一做 BH 校正、样本量过滤有效样本少于 30 的对直接跳过、以及把结果写成带方法名和置信区间的表方便回溯。还有一点经验之谈自动化分发的结果一定要抽查。我做过一次自动跑完 200 个变量对发现有几十对的 p 值全是 1.0查下来是一列常量某个渠道在某个时间段没有数据方差为 0检验直接失效。没有异常检测的自动化流程出错的时候是静默的。7. 常见问题速查与排查思路方法学完了剩下的就是踩坑。这一节把我在实际项目里反复遇到的问题整理成表附上判断依据和处置方式。7.1 问题现象与对应处置表现象可能原因判断依据处置方式Pearson 显著但 Spearman 不显著异常值或非线性散点图、箱线图查异常值或改报 Spearman相关系数极小但 p 极小样本量过大看 n 和置信区间报效应量与区间别只报 p卡方显著但看不出差异期望频数或样本量问题检查交叉表百分比报 Cramérs V 和百分比分布ANOVA 显著但事后比较全不显著未校正 / 组数太多看比较次数用 BH 校正或减少比较两组均值差很大但 t 检验不显著组内方差极大 / 样本小看标准差和 n增大样本或改用非参数点二列相关系数偏低两组样本量悬殊看两组比例改用 Cohens d 判断二分类编码后相关方向反了编码顺序与业务方向不一致看 0/1 对应什么明确编码约定并写入文档多个变量两两比较出现大量显著多重比较膨胀比较次数FDR 校正或降维后再看表里第三条我要多嘴一句。卡方显著但事后看不出差异八成是交叉表的百分比差异很小、但样本量大到能检出。这种时候别硬讲有关系直接写存在统计关联但效应量很小V 0.04业务上可视为无实质差异比编一个故事强。第七条也很常见。is_member如果 1 表示非会员算出来的相关系数符号就是反的解释的时候会得出完全相反的结论。这个错误极其低级但极其常见尤其是多人协作、字段含义靠口口相传的项目。我的做法是在数据字典里强制写取值含义这一列而且要求写具体1 是会员不是1 是。7.2 我踩过的几个坑第一个坑是把相关当因果用。这个错误不是统计层面的是表达层面的。有一次做渠道与客单价的分析算出来推荐渠道的客单价显著更高报告里写了一句推荐渠道带来更高客单价被业务方理解成多投推荐渠道就能提升客单价。实际上很可能是高价值用户更倾向于走推荐路径方向完全相反。后来我改了个规矩所有相关性结论只写关联伴随差异动词一律不用导致带来提升。这个习惯救过我好几次。第二个坑是忽略辛普森悖论。分组看和整体看结论相反的情况在分类对连续的分析里特别容易发生。比如整体看会员的客单价更高但按城市分层之后每个城市里会员的客单价都更低——原因可能是会员在高线城市占比更高而高线城市本身客单价就高。这个坑的排查方式很朴素算完整体一定要按关键维度分层再算一遍两层结论不一致就说明有混杂变量得用分层分析或者回归控制。第三个坑是样本量不均衡导致的误判。多分类对连续一组 400 人一组 20 人ANOVA 的 F 值主要被大组驱动小组的均值极不稳定却参与了组间平方和的计算。这种情况我一般会做两件事一是看每组的置信区间小组的区间宽得离谱就直接标注样本不足结论仅供参考二是补一个 bootstrap 的组间差异估计看看重抽样下结论稳不稳。第四个坑是测量工具的刻度问题。问卷里的李克特量表从非常不满意到非常满意标 1 到 5看起来是等距的但受访者对满意和非常满意的心理距离跟一般和满意的距离往往不一样。这种不等距在算 Pearson 时会引入偏差算 Spearman 就没问题。所以但凡数据来源是主观量表我现在默认先跑 SpearmanPearson 只作为对照。第五个坑是把统计不显著当成没有关系。这个在样本量小的时候特别危险。n 40 的时候即使真实相关有 0.3检出的概率也不到一半。所以不显著时正确的表述是在当前样本量下未检出显著关联而不是两者无关。如果业务上这个问题很重要那就该去算需要多大样本量才能检出预期效应而不是拿一个不显著的结果草草收场。这套东西用熟了之后回头看那些直接df.corr()一把梭的做法会觉得挺粗糙的。但说实话大多数业务场景也不需要多精细——知道该看哪个指标、知道 p 值不能单独用、知道分组之后要再检查一遍就已经能避开 80% 的坑了。剩下的 20%等你真的被一个辛普森悖论坑过一次自然就记住了。
应用指南)
