回归斜率的非参数统计显著性测试问题求助

问题分析 #

你当前遇到的核心问题有两个：

1. scipy.stats.mannwhitneyu的默认行为导致输出维度不匹配：该函数会将多维输入按最后一维做批量检验，因此输入(721,1440)的两个数组时，会对每一列的721个元素做检验，最终输出1440个结果（一维数组），而非每个网格点单独的检验值。
2. 方法选择偏差：Mann-Whitney U检验的核心是对比两个独立样本的分布差异，并不适用于单个回归斜率与0的显著性检验。若要检验每个网格点的时间序列趋势（回归斜率）是否显著，Mann-Kendall趋势检验才是气象/海洋领域常用的非参数方法，完全符合审稿人要求的“非参数、可与Student's T-test对比”的条件。

解决方案 #

方案1：使用Mann-Kendall检验（推荐） #

Mann-Kendall检验直接针对时间序列的单调趋势做显著性检验，无需假设正态分布，可直接输出每个网格点的p值，维度与原始斜率数组一致。

使用scipy.stats.mannkendall结合向量化操作实现：

import numpy as np
from scipy.stats import mann_kendall

# 假设原始时间序列数据是shape=(n_time, 721, 1440)的数组data_y（大气变量）
# 初始化p值数组
p_values = np.zeros((721, 1440))

# 遍历每个网格点计算检验结果
for i in range(721):
    for j in range(1440):
        _, p = mann_kendall(data_y[:, i, j])
        p_values[i, j] = p

# 若需同时计算非参数斜率（替代线性回归斜率，可选）
from scipy.stats import theilslopes
slopes = np.zeros((721, 1440))
for i in range(721):
    for j in range(1440):
        slope, _, _, _ = theilslopes(data_y[:, i, j])
        slopes[i, j] = slope

方案2：修正Mann-Whitney U检验的使用逻辑（若严格遵循审稿人建议） #

若必须使用Mann-Whitney U检验，需明确检验逻辑：为每个网格点构造“零斜率的模拟样本”（如基于随机打乱的时间序列计算的斜率），再与实际斜率做分布对比。以下是针对每个网格点单独计算的实现：

import numpy as np
from scipy.stats import mannwhitneyu

# 假设实际回归斜率数组regression的shape=(721,1440)
# 为每个网格点构造零分布样本（示例：生成1000个零附近的随机样本，模拟无趋势情况）
n_sim = 1000
zero_samples = np.random.normal(loc=0, scale=np.std(regression), size=(n_sim, 721, 1440))

# 初始化p值数组
p_values = np.zeros((721, 1440))

# 遍历每个网格点做检验
for i in range(721):
    for j in range(1440):
        # 对比实际斜率与零分布样本
        _, p = mannwhitneyu([regression[i,j]]*n_sim, zero_samples[:,i,j], alternative='two-sided')
        p_values[i,j] = p

快速向量化优化 #

如果嵌套循环速度较慢，可使用numpy.apply_along_axis实现向量化加速：

def mk_test(arr):
    _, p = mann_kendall(arr)
    return p

# data_y shape=(n_time, 721, 1440)
p_values = np.apply_along_axis(mk_test, 0, data_y)
# 输出p_values的shape=(721,1440)，与需求完全匹配

内容的提问来源于stack exchange，提问作者Megan Franke