如何对时间序列遥感数据应用CDF匹配法进行归一化处理？

核心需求说明 #

需要将卫星遥感土壤湿度（SSM）数据通过CDF匹配法归一化，使其时间序列的均值、方差、偏度、峰度与地面观测参考数据（P）近似一致，为后续可靠性验证做准备。

现有代码问题分析 #

1. 自定义prob_data函数未关联排序后的数据，无法正确构建完整的CDF曲线
2. 插值逻辑偏离CDF匹配核心：标准CDF匹配是将有偏数据的累积概率映射到参考数据的对应分位值，而非反向插值
3. 额外的多项式拟合属于冗余步骤，会引入不必要的拟合误差
4. 绘图x轴标注错误（应为变量值而非时间）
5. 数据处理顺序不合理：先归一化再删除空值，会导致归一化结果受空值干扰
6. 最终数据合并时存在索引不匹配风险

修正后的完整代码 #

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy.interpolate import interp1d
import scipy.stats

# 1. 数据读取与预处理
data = pd.read_excel('Data.xlsx', sheet_name='l4')
data['Timestamp'] = pd.to_datetime(data['Timestamp'])
data.set_index('Timestamp', inplace=True)

# 先删除空值，避免后续计算受无效值影响
data = data.dropna(subset=['P', 'SSM'])

# 可选归一化函数（若原始数据量纲差异大时使用）
def normalize_data(value):
    min_val = np.min(value)
    max_val = np.max(value)
    return (value - min_val) / (max_val - min_val)

# 提取参考数据（地面观测P）和有偏数据（卫星SSM）
ref_data = data['P'].values
biased_data = data['SSM'].values

# 可选：若需要先做0-1归一化，取消注释以下两行
# ref_data = normalize_data(ref_data)
# biased_data = normalize_data(biased_data)

# 2. 计算经验CDF（累积分布函数）
def compute_cdf(data):
    sorted_data = np.sort(data)
    n = len(sorted_data)
    # 采用(n)/(n+1)计算累积概率，避免概率值达到1
    cdf_prob = np.arange(1, n+1) / (n + 1)
    return sorted_data, cdf_prob

# 生成参考数据和有偏数据的CDF
sorted_ref, cdf_ref = compute_cdf(ref_data)
sorted_biased, cdf_biased = compute_cdf(biased_data)

# 3. CDF匹配核心：构建概率-值映射关系
# 用参考数据的CDF构建插值函数，实现概率到参考值的映射
interp_func = interp1d(cdf_ref, sorted_ref, bounds_error=False, fill_value=(sorted_ref[0], sorted_ref[-1]))

# 对单个数据点执行CDF匹配：找到其在自身CDF中的概率，再映射到参考数据的对应值
def cdf_match_single(x, sorted_data, cdf_prob):
    idx = np.searchsorted(sorted_data, x, side='right')
    if idx == 0:
        prob = cdf_prob[0]
    elif idx == len(sorted_data):
        prob = cdf_prob[-1]
    else:
        # 线性插值计算当前值对应的累积概率
        prob = cdf_prob[idx-1] + (x - sorted_data[idx-1])*(cdf_prob[idx]-cdf_prob[idx-1])/(sorted_data[idx]-sorted_data[idx-1])
    return interp_func(prob)

# 批量处理所有卫星数据
corrected_data = np.array([cdf_match_single(x, sorted_biased, cdf_biased) for x in biased_data])

# 4. 统计特征验证（检查匹配效果）
def get_stats(data):
    return {
        '均值': round(np.mean(data), 4),
        '方差': round(np.var(data), 4),
        '偏度': round(scipy.stats.skew(data), 4),
        '峰度': round(scipy.stats.kurtosis(data), 4)
    }

print("参考数据统计特征：", get_stats(ref_data))
print("原始卫星数据统计特征：", get_stats(biased_data))
print("匹配后卫星数据统计特征：", get_stats(corrected_data))

# 5. 绘图展示CDF匹配效果
plt.figure(figsize=(8,6))
plt.plot(sorted_ref, cdf_ref, 'k-', linewidth=1.5, label='参考观测数据')
plt.plot(sorted_biased, cdf_biased, 'k--', linewidth=1.5, label='原始卫星数据')
plt.plot(np.sort(corrected_data), compute_cdf(corrected_data)[1], 'k:', linewidth=1.5, label='匹配后卫星数据')

plt.xlabel('土壤湿度值')
plt.ylabel('累积分布概率（CDF）')
plt.legend(loc='best')
plt.gcf().set_facecolor('white')
plt.savefig('CDF_SSM_Match.png', format='png', dpi=300)
plt.show()

# 6. 结果保存
data_correct = pd.DataFrame({
    'Timestamp': data.index,
    'P': ref_data,
    'SSM': biased_data,
    'SSM_cdf_matched': corrected_data
})

# 可选：重新索引到每日时间序列（补全缺失日期）
data_correct = data_correct.set_index('Timestamp').reindex(
    pd.date_range('2016-01-01', '2021-12-31', freq='D')
).rename_axis(['Date']).reset_index()

data_correct.to_excel('L4_CDF_Matched_Data.xlsx', index=False)

关键步骤说明 #

1. 数据预处理：优先删除空值，避免无效数据干扰后续计算；归一化步骤可根据原始数据量纲差异选择是否启用。
2. CDF计算：对排序后的数据计算经验累积概率，确保每个数据点对应唯一的概率位。
3. CDF映射：核心逻辑是将卫星数据的每个值映射到参考数据的同概率分位值，确保两者分布特征完全匹配。
4. 效果验证：通过统计指标直接验证匹配效果，确认均值、方差等特征已与参考数据对齐。
5. 结果保存：按时间索引整理数据，方便后续可靠性验证分析。

内容的提问来源于stack exchange，提问作者Sofia Ortenzi