水文分析与计算——年均流量趋势检验(Mann-Kendall法、线性回归法)

[cpp]
//年均流量趋势检验.h

//年均流量Mann-Kendall法趋势分析
void MannKendall()
{
    using namespace std;
    int S = 0;//检验的统计变量
    double VarS,//统计变量S的方差
        Z ;//标准正态统计变量方差
    S = 0;
    for(int i = 0; i < Y; i++)
        for(int j = i + 1; j < Y; j++)
        {
            if(YearQ[j]>YearQ[i]) S++;
            if(YearQ[j]<YearQ[i]) S--;
        }
    VarS = 0;
    VarS = Y*(Y - 1)*(2*Y + 5)/18.0;
    if(S > 0) Z = (S - 1)/pow(VarS, 0.5);
    if(S < 0) Z = (S + 1)/pow(VarS, 0.5);
    cout<<"年均流量趋势检验——Mann-Kendall检验："<<endl
        <<"标准正态统计变量Z= "<<Z<<endl
        <<"Mann-Kendall检验通过请输入1,否则请关闭！"<<endl;
    cin>>Z;//控制台停留
    cout<<endl;
}

//年均流量线性回归法法趋势分析
double Normal(double z)
{//返回标准正态分布的密度函数
    double temp;
    temp=exp((-1)*z*z/2)/sqrt(2*PI);
    return temp;
}
double NormSDist(const double z)
{//返回标准正态分布的累积频率函数
    if(z > 6) return 1;
    if(z < -6) return 0;
    static const double gamma = 0.231641900,
        a1 = 0.319381530,
        a2 = -0.356563782,
        a3 = 1.781477973,
        a4 = -1.821255978,
        a5 = 1.330274429;
    double k = 1.0 / (1 + fabs(z) * gamma);
    double n = k * (a1 + k * (a2 + k * (a3 + k * (a4 + k * a5))));
    n = 1 - Normal(z) * n;
    if(z < 0)
        return 1.0 - n;
    return n;
}
void XianXingJianYan()
{//线性回归检验
    using namespace std;
    double AverageYearQ = 0,
        AverageT = (1+Y)/2.0,
        a, b,//待定回归系数
        r = 0,//线性相关系数
        t,//t统计量
        sigmaT =0,
        sigmaYearQ = 0,//均方差
        NewYearQ[Y],//按升序排列的年均流量
        Fn,//样本累积频率
        F0,//理论累积频率
        D_n_alpha = 0.202737,//显著水平为alpha且样本容量为n时的拒绝临界值
        MaxD = 0,//max(|Fn - F0|)
        temp,
        sigmab = 0;//回归系数b标准方差
// int order;//升序排序年均流量
    for(int i = 0; i < Y; i++)
    {
        AverageYearQ += YearQ[i];
    }
    AverageYearQ /= Y;
    for(int i = 0; i < Y; i++)
    {
        r += (i - AverageT)*(YearQ[i] - AverageYearQ);
        sigmaT += pow(i - AverageT, 2);
        sigmaYearQ += pow(YearQ[i] - AverageYearQ, 2);
    }
    r /= pow(sigmaT*sigmaYearQ, 0.5);
    sigmaT = pow(sigmaT/(Y - 1), 0.5);
    sigmaYearQ = pow(sigmaYearQ/(Y - 1), 0.5);
    for(int i = 0; i < Y; i++)
        NewYearQ[i] = YearQ[i];
    for(int i = 0; i < Y - 1; i++)
    {
        for(int j = i + 1; j < Y; j++)
            if(NewYearQ[i] > NewYearQ[j])
            {
                temp = NewYearQ[i];
                NewYearQ[i] = NewYearQ[j];
                NewYearQ[j] = temp;
            }
    }
    for(int i = 0; i < Y; i++)
    {
        Fn = (double)(i+1)/(Y + 1);
        F0 = NormSDist((NewYearQ[i] - AverageYearQ)/sigmaYearQ);
        if(MaxD < fabs(Fn - F0)) MaxD = fabs(Fn - F0);
    }
    cout<<"年均流量趋势检验——线性回归检验："<<endl
        <<"正态分布K-S检验统计量D ="<<MaxD<<endl
        <<"K-S检验拒绝临界值D(n, a)="<<D_n_alpha<<endl;
    b = r*sigmaYearQ/sigmaT;
    a = AverageYearQ - b*AverageT;
    for(int i = 0; i < Y; i++)
        sigmab += pow(YearQ[i] - (a + b*i), 2);
    sigmab = pow(sigmab/(Y - 2), 0.5)/(pow(sigmaT, 2)*(Y - 1));
    t = b/sigmab;
    cout<<"线性相关系数r = "<<r<<endl
        <<"年均流量Q倚时序t的回归系数估计值分别为："<<endl
        <<"a = "<<a<<endl
        <<"b = "<<b<<endl
        <<"假设检验统计量t = "<<t&l

补充：软件开发 , C++ ,