气象统计 第三次实习 落后交叉相关系数和偏相关系数

气象统计 第三次实习

一、实习题目

二、实习程序

import pandas as pd
import numpy as np

##读文件
data1=np.array(pd.read_table('mon1.txt',header=None))
data2=np.array(pd.read_table('mon2.txt',header=None))
data12=np.array(pd.read_table('mon12.txt',header=None))
##简单相关系数,注意x,y长度相等
def xgxs(x,y):
    n1=len(x)    
    ave1=np.mean(x)
    ave2=np.mean(y)
    s1= np.std(x)
    s2= np.std(y)
    sum=0
    for i in range(0,n1):
        sum=(x[i]-ave1)*(y[i]-ave2)+sum
    r=sum/s1/s2/n1
    return r 
    
###落后交叉相关系数函数,注意x,y长度相等
def lhjcxgxs(x,y):
    n1=len(x)
    N=int(n1/2)
    ave1=np.mean(x)
    ave2=np.mean(y)
    s1= np.std(x)
    s2= np.std(y)
    r=np.zeros(N)
    for i in range(0,N):
        sum=0
        for j in range(0,n1-i-1):
            sum=((x[j]-ave1)/s1)*((y[i+1+j]-ave2)/s2)+sum
        r[i]=1/(n1-i-1)*sum
    return r

##12月与1月的落后交叉相关系数
print('12月与1月的落后交叉相关系数')
r_xy=lhjcxgxs(data12,data1)
for i in range(0,10):
    print('t=%d  r_xy=%lf'%((i+1),r_xy[i]))
##12月与1月的偏相关系数
r=(xgxs(data12,data1)-xgxs(data1,data2)*xgxs(data12,data2))/((1-xgxs(data1,data2)**2)*(1-xgxs(data12,data2)**2))**0.5
print('12月与1月的偏相关系数:%lf'%(r))
##12月与2月的落后交叉相关系数
print('12月与2月的落后交叉相关系数')
r_xy=lhjcxgxs(data12,data2)
for i in range(0,10):
    print('t=%d  r_xy=%lf'%((i+1),r_xy[i]))
##12月与2月的偏相关系数
r=(xgxs(data12,data2)-xgxs(data1,data2)*xgxs(data12,data1))/((1-xgxs(data12,data1)**2)*(1-xgxs(data1,data2)**2))**0.5
print('12月与2月的偏相关系数:%lf'%(r))

三、例图

四、代码及源文件

实习三

程序中如有疏漏欢迎指正