气象学中的合成分析

一、前言

在复杂多变的气象系统中，合成分析（Composite Analysis） 是一种揭示天气事件背后共性规律的核心研究方法。它通过将大量相似天气或气候事件的数据进行空间-时间的整合与对比，抽离出关键信号，从而帮助科学家理解极端天气成因、气候模态的演变机制，甚至改进数值预报模型的精度。

二、原理

合成分析（Composite Analysis）是气象学中揭示天气或气候事件共性规律的核心统计工具，其本质是通过对同类事件的时空数据叠加与对比，提取共性信号并探索物理机制。

• 事件定义：基于物理标准（如台风强度、热浪阈值）或统计阈值（如降水距平标准差），筛选出满足条件的同类事件样本（如50次极端降雨事件）。
• 数据对齐：以事件发生的时间或空间关键点（如台风登陆时刻、阻塞高压中心）为基准，对齐所有样本的时空数据场（如气压、温度、风场）。
• 信号提取：通过计算均值场或异常场（合成场），消除随机波动（噪音），凸显事件共性特征。

三、函数

import numpy as np
from datetime import datetime
from scipy.stats.mstats import ttest_ind
def composite_analysis(data: np.ndarray,  index: np.ndarray,  date: np.ndarray,  threshold: float = 0.8) -> tuple:
    '''
    气象合成分析函数（修正版）
    
    参数说明：
    data   - 三维气象场数据，形状为(time, lat, lon)
             type: numpy.ndarray 
             example: 温度场数据，维度为[365, 180, 360]表示365天x180纬度x360经度
    index  - 一维气候指数序列，形状为(time,)
             type: numpy.ndarray 
             example: 标准化后的ENSO指数序列
    date   - 一维日期序列，形状为(time,)
             type: numpy.ndarray 
             example: numpy.datetime64 类型数组，如['2000-01-01', ...]
    threshold - 阈值系数，默认0.8个标准差
             type: float
             range: 建议0.5-1.5之间

    返回：
    composite_high - 高指数合成场异常（相对于气候态），形状(lat, lon)
    composite_low  - 低指数合成场异常（相对于气候态），形状(lat, lon)
    date_high      - 高指数事件日期数组
    date_low       - 低指数事件日期数组

    算法原理：
    1. 对指数进行标准化处理（Z-score）
    2. 筛选超过阈值的高/低指数事件日
    3. 计算事件日的平均气象场
    4. 计算相对于气候平均态的异常场
    '''
    
    # 输入验证（新增）
    if data.shape[0]  != index.size  or data.shape[0]  != date.size: 
        raise ValueError("输入维度不匹配：data.shape[0]={},  index.size={},  date.size={}" 
                        .format(data.shape[0],  index.size,  date.size)) 
    if not isinstance(threshold, (int, float)):
        raise TypeError("阈值参数应为数值类型")
    
    # 指数标准化（修正为样本标准差）
    index_mean = np.nanmean(index)   # 忽略NaN值计算
    index_std = np.nanstd(index,  ddof=1)  # 使用样本标准差（ddof=1）
    index_stdized = (index - index_mean) / index_std
    
    # 事件筛选（新增边界条件检查）
    high_mask = index_stdized > threshold
    low_mask = index_stdized <= -1*threshold
    
    # 异常处理：检查有效事件数量（新增）
    if np.sum(high_mask)  < 3:  # 至少需要3个样本保证统计显著性
        raise ValueError(f"高指数事件样本不足（{np.sum(high_mask)} 个），请降低阈值")
    if np.sum(low_mask)  < 3:
        raise ValueError(f"低指数事件样本不足（{np.sum(low_mask)} 个），请升高阈值")

    # 数据提取
    date_high = date[high_mask]
    date_low = date[low_mask]
    data_high = data[high_mask]
    data_low = data[low_mask]
    
    # 气候态计算（优化NaN处理）
    clim_mean = np.nanmean(data,  axis=0)  # 沿时间维度的气候平均场
    # _, p_high = ttest_ind(data_high, clim_mean, equal_var=False)
    # _, p_low = ttest_ind(data_low, clim_mean, equal_var=False)
    # # 原错误代码段
    # _, p_high = ttest_ind(data_high, clim_mean, equal_var=False)
    # _, p_low = ttest_ind(data_low, clim_mean, equal_var=False)

    # 修正后代码
    from scipy.stats  import ttest_1samp  # 注意导入路径变化

    _, p_high = ttest_1samp(data_high, popmean=clim_mean, axis=0)
    _, p_low = ttest_1samp(data_low, popmean=clim_mean, axis=0)
    print(data_high.shape)
    print(data_low.shape)
    # 合成分析计算（新增加权平均选项）
    composite_high = np.nanmean(data_high,  axis=0) - clim_mean
    composite_low = np.nanmean(data_low,  axis=0) - clim_mean
    
    return composite_high, composite_low, p_high, p_low, date_high, data_low

四、参考文献

1. https://www.heywhale.com/mw/project/61e4e987ddda3c0017b4b1d4