绝对准确的计算STI的Python代码

___Y1

已于 2024-05-23 12:37:03 修改

阅读量359

点赞数 5

分类专栏：气象数据处理文章标签： python numpy 机器学习

于 2024-05-15 17:56:23 首次发布

本文链接：https://blog.csdn.net/weixin_44259522/article/details/138917631

版权

气象数据处理专栏收录该内容

47 篇文章 0 订阅

订阅专栏

这个脚本看起来是用来处理包含气象数据的Excel文件，并计算相关的统计量。让我来解释每个函数的作用：

1. `calculate_alpha_beta(x)`:
- 这个函数计算了参数alpha和beta，这是用来拟合Gamma分布的参数。
- 参数x是一个包含气象数据的数组。
- 这个函数首先计算了x的均值x_bar，然后根据公式计算了alpha_hat和beta_hat。

举例：假设输入的x是一个包含某地区一年的温度数据的数组，这个函数会根据这些数据计算出适合描述这个温度分布的Gamma分布的参数alpha和beta。

2. `calculate_H(x, alpha, beta)`:
- 这个函数用来计算累积分布函数（CDF），即用参数alpha和beta拟合的Gamma分布在给定的数据点x处的值。
- 参数x是一个包含气象数据的数组，而alpha和beta是Gamma分布的参数。
- 这个函数使用了SciPy库中的gamma.cdf函数来计算Gamma分布的CDF值。

举例：如果我们知道某地区一年的温度数据符合Gamma分布，我们可以使用这个函数来计算在某个特定温度值处的累积分布函数值。

3. `standardize_H(G_x)`:
- 这个函数用来标准化给定的累积分布函数值。
- 参数G_x是累积分布函数的值。
- 这个函数将累积分布函数值转换为标准正态分布的z-score值。

举例：如果我们有某地区一年的温度数据的累积分布函数值，我们可以使用这个函数将这些值标准化为标准正态分布的z-score值，以便进行比较或其他统计分析。

4. `process_file(file_path)`:
- 这个函数是整个处理过程的主函数，它读取Excel文件中的数据，然后依次调用前面的函数来计算统计量并保存结果。
- 参数file_path是要处理的Excel文件的路径。

举例：假设我们有一个Excel文件包含了某地区多年的温度数据，这个函数会读取这个文件，然后计算出每年的STI（标准温度指数），最后将结果保存到一个新的Excel文件中。

import numpy as np
from scipy.stats import gamma
import pandas as pd
import os

def calculate_alpha_beta(x):
    n = len(x)
    x_bar = np.mean(x)
    A = np.log(x_bar) - (1/n) * np.sum(np.log(x))
    alpha_hat = (1 + np.sqrt(1 + 4 * A / 3)) / (4 * A)
    beta_hat = x_bar / alpha_hat
    return alpha_hat, beta_hat

def calculate_H(x, alpha, beta):
    G_x = gamma.cdf(x, alpha, scale=beta)
    return G_x

def standardize_H(G_x):
    z_score = (G_x - np.mean(G_x)) / np.std(G_x)
    STI = z_score
    return STI

def process_file(file_path):
    data = pd.read_excel(file_path)
    x = data['TMED'].values
    alpha_hat, beta_hat = calculate_alpha_beta(x)
    G_x = calculate_H(x, alpha_hat, beta_hat)
    STI = standardize_H(G_x)
    data['STI'] = STI
    output_file_path = os.path.splitext(file_path)[0] + '_sti.xlsx'
    data[['year', 'STI']].to_excel(output_file_path, index=False)
    print("STI结果已保存到:", output_file_path)

folder_path = "D:\copula\新疆年气象数据"

for file_name in os.listdir(folder_path):
    if file_name.endswith('.xlsx'):
        file_path = os.path.join(folder_path, file_name)
        process_file(file_path)

这个脚本涉及到了一些统计量的计算，主要涉及到Gamma分布和标准化的过程。下面是相关的计算公式总结：

1. **计算Alpha和Beta的公式**：
- Alpha和Beta是用来拟合Gamma分布的参数。
- Alpha的计算公式为：
$\hat{\alpha} = \frac{{1 + \sqrt{1 + \frac{4A}{3}}}}{{4A}}$
- Beta的计算公式为：
$\hat{\beta} = \frac{{\bar{x}}}{{\hat{\alpha}}}$
其中：
- $A = \log(\bar{x}) - \frac{1}{n} \sum \log(x)$ ，$\bar{x}是样本均值，\(n$ 是样本大小，$x$ 是样本数据。