最近在看中长期水文预报,打算使用python语言实现课本的模型并进行实例的计算结果的检验,为了监督自己和整理记录自己实现的模型代码,打算写博客记录自己的代码实现和部分思路。首先,自己不使用现成的模块实现的是一元线性回归模型,然后和模块实现对比,学习模块的调用。接下来是直接编写的程序文件,并计算了几个检验系数,没有使用sklearn模块:
#导模块
import numpy as np
import pandas as pd
#全局变量初始化
X=[]
Y=[]
Y_predict=[]
S_xy=0
S_xx=0
S_yy=0
U=0
#从csv文件中读取数据
def get_data(file_name):
#1.利用pandas读取csv
data=pd.read_csv(file_name,header=None)
#2.赋值预报因子和预报对象并输出
X=data[0]#预报因子
Y=data[1]#预报对象
return X,Y
#print(get_data('ycfj.csv'))
X,Y=get_data('ycfj.csv')
#计算预报因子和预报对象的均值
X_mean=np.mean(X)
Y_mean=np.mean(Y)
#计算回归系数
for i in range(len(X)):
S_xy+=(X[i]-X_mean)*(Y[i]-Y_mean)
S_xx+=pow(X[i]-X_mean,2)
S_yy+=pow(Y[i]-Y_mean,2)
b1=S_xy/S_xx
b0=Y_mean-b1*X_mean
print('b1:'+str(b1)+' b0:'+str(b0))
#t检验(显著性检验)
t=S_xy*pow((len(X)-2)/(S_xx*S_yy-pow(S_xy,2)),0.5)
print(t)
#相关系数检验
r=S_xy/pow(S_xx*S_yy,0.5)
#可决系数检验
for i in range(len(X)):
y_predict=b0+b1*X[i]
Y_predict.append(y_predict)
U+=pow(Y_predict[i]-Y_mean,2)
R2=U/S_yy
print(R2)
下面的是使用了sklearn模块,并进行了绘图:
#导模块
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.linear_model import LinearRegression
#全局变量初始化
X=[]
Y=[]
Y_predict=[]
S_xy=0
S_xx=0
S_yy=0
U=0
#从csv文件中读取数据
def get_data(file_name):
#1.利用pandas读取csv
data=pd.read_csv(file_name,header=None)
#2.赋值预报因子和预报对象并输出
X=data[0]#预报因子
Y=data[1]#预报对象
return X,Y
#print(get_data('ycfj.csv'))
X,Y=get_data('ycfj.csv')
X=np.array(X).reshape(-1,1)#将数组转换?
Y=np.array(Y).reshape(-1,1)
def my_linear_model(X,Y):
#1.调用模块,构造回归对象
regr=LinearRegression()
regr.fit(X,Y)
#2.构造返回字典
predictions={}
#2.1截距
predictions['intercept']=regr.intercept_
#2.2斜率值
predictions['coefficient']=regr.coef_
#2.3预测值
predictions['predict_value']=regr.predict(X)
return predictions
# 绘图
def show_linear_line(X,Y):
predictions=my_linear_model(X,Y)
# 1.绘出已知数据散点图
plt.scatter(X,Y,color='blue')
# 2. 绘出预测直线
plt.plot(X,predictions['predict_value'],color ='red',linewidth=4)
plt.title('Monolinear Regression Analysis')
plt.xlabel('X')
plt.ylabel('Y')
plt.show()
def main():
show_linear_line(X,Y)
if __name__=='__main__':
main()