水文分析与计算——年均流量趋势检验(Mann-Kendall法、线性回归法)
//年均流量趋势检验.h//年均流量Mann-Kendall法趋势分析void MannKendall(){using namespace std;int S = 0;//检验的统计变量doubleVarS,//统计变量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<<endl<<"线性回归检验通过请输入1,否则请关闭！"<<endl;cin>>t;//控制台停留cout<<endl<<endl<<endl;}