Python实现曼-肯德尔检验法

直孔吧-坚扎

'''
实现曼-肯德尔检验法encoding=UTF-8'''
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from matplotlib.ticker import MultipleLocator, FormatStrFormatter

## 数据处理部分
def U(tm):
    ## 计算UF
    n=len(tm)
    s1= np.arange(0, 91, 1)
    r=np.zeros(shape=(91,91))
    for i in range(n):
        for j in range(n):
            if tm>tm[j]:
                r[j]=1
            else:
                r[j]=0
    s1[0]=0.0
    for i in range(n-1):
        i+=1
        s1=s1[i-1]
        for j in range(i):
            s1+=r[j]
    E1 = []
    var1 = []
    for k in range(n):
        k += 1.0
        jz = (k * (k - 1.0)) / 4.0
        fc = (k * (k - 1.0) * (2.0 * k + 5.0)) / 72.0
        E1.append(jz)
        var1.append(fc)

    UF = []
    for k in range(n):
        UF.append((s1[k] - E1[k]) / np.sqrt(var1[k]))
    UF[0]=0.0

    ## 计算UB
    s2 = np.arange(0, 91, 1)
    s2[n-1] = 0.0
    i=n-2
    while (i>=0):   
        s2=s2[i+1]
        j=n-1
        while (j>i):
            s2+=r[j]
            j -= 1
        i-=1

    E2 = []
    var2 = []
    k=n
    while (k>=0):
        jz = (k * (k - 1.0)) / 4.0
        fc = (k * (k - 1.0) * (2.0 * k + 5.0)) / 72.0
        k-=1
        E2.append(jz)
        var2.append(fc)

    ## 计算UB
    UB = []
    i=0
    while (i<n):
        UB.append((E2-s2) / np.sqrt(var2))
        i += 1
    UB[n-1]=0.0
    UF=pd.DataFrame(UF)
    UB=pd.DataFrame(UB)
    return UF,UB

## 画图部分
def plot(year,uf,ub):

    year=np.array(year)
    uf=np.array(uf)
    ub=np.array(ub)
    y=np.arange(-2,9)

    fig,axs=plt.subplots(1,figsize=(9,4))
    axs.set_title('Mann-Kendall')
    axs.set_xlabel('year')
    axs.set_ylabel('statistical magnitude')
    xmajorLocator = MultipleLocator(10);
    xmajorFormatter = FormatStrFormatter('%1d')
    xminorLocator = MultipleLocator(10)

    ymajorLocator = MultipleLocator(1)
    ymajorFormatter = FormatStrFormatter('%1.1f')
    yminorLocator = MultipleLocator(8)

    axs.xaxis.set_major_locator(xmajorLocator)
    axs.xaxis.set_major_formatter(xmajorFormatter)
    axs.yaxis.set_major_locator(ymajorLocator)
    # axs.yaxis.set_major_formatter(ymajorFormatter)

    axs.xaxis.set_minor_locator(xminorLocator)
    axs.yaxis.set_minor_locator(yminorLocator)

    plt.ylim(-2,6)
    axs.plot(year,uf, c='b', ls="-", linewidth=1.5)
    axs.plot(year, ub, c='b', ls="--", linewidth=1.5)
    axs.legend(["UF", "UB", loc="upper right")

    plt.savefig('./data/Kendall.png')
    # plt.show()

def main():
    path='./data/ysj_data.xlsx'
    data=pd.read_excel(path,skiprows=0,names=['YEAR','TM'])
    year_data=data['YEAR']
    tm_data=data['TM']
    uf_data,ub_data = U(tm_data)
    plot(year_data,uf_data,ub_data)

if __name__ == '__main__':
    main();

lancerv8

第一次看到Python做的MK检验程序，学习一下。

rgb001

感谢分享！

kmdqlb

老师，程序所用数据xlsx也请传上来，让大家都学习一下

勿忘草

学习了，请问能不能将数据文件也放上来