随便记录一下,财政收入的数据搜一下就有了。
本文主要是三种回归的对比:岭回归、lasso回归、弹性网络回归。
ps:不要问问题,因为我也不知道答案(如果有人问的话)
描述性分析
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
from sklearn.linear_model import LinearRegression
import statsmodels.api as sm
sns.set() #if you want to use seaborn themes with matplotlib functions
import warnings
warnings.filterwarnings('ignore')
rand_state= 1000
df=pd.read_csv(r'D:\桌面\python代码+数据\chapter6\demo\data\data.csv')
df.head()
sns.scatterplot(x='x1', y='y', data=df)
plt.show()
#在训练之前,先进行标准化
from sklearn.preprocessing import StandardScaler
scaler=StandardScaler()
df_sc= scaler.fit_transform(df)
df_sc[0]
df.describe()
sns.scatterplot(x='x2', y='y', data=df_sc)
plt.show()
划分数据
y = df_sc['y']
X = df_sc.drop('y', axis=1) # becareful inplace= False
from sklearn.model_selection import train_test_split
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=rand_state)
X_train.head()
从上表中可以看出,变量彼此高度相关。 让我们首先运行线性回归(使用 statsmodels.api),并将结果视为基准。
X_test_wc = sm.add_constant(X_test)
X_train_wc = sm.add_constant(X_train)
热力图
corrdf=df.corr()
corrdf
plt.figure(figsize=(10,10))
sns.heatmap(corrdf,annot=True,cmap='RdGy')
plt.show()
label=['y']
gdp=(list(df.y))
plt.figure(figsize=(6,4))
plt.boxplot(gdp,notch=True,labels=label,meanline=True)
plt.show()
OLS回归的结果
model = sm.OLS(y_train, X_train_wc).fit()#添加截距列
model.summary()
from sklearn.linear_model import LinearRegression, RidgeCV, Lasso, LassoCV, ElasticNet, ElasticNetCV
from sklearn.linear_model import Ridge
model_linear = LinearRegression()
model_ridge = Ridge(alpha=10)
model_lasso = Lasso(alpha=10)
model_net = ElasticNet(alpha=10)
y_hat_linear= model_linear.fit(X_train, y_train).predict(X_test)
y_hat_ridge = model_ridge.fit(X_train, y_train).predict(X_test)
y_hat_lasso = model_lasso.fit(X_train, y_train).predict(X_test)
y_hat_net = model_net.fit(X_train, y_train).predict(X_test)
df_predictions = pd.DataFrame({'y_test':y_test,
'y_hat_linear':y_hat_linear,
'y_hat_ridge':y_hat_ridge,
'y_hat_lasso':y_hat_lasso,
'y_hat_net':y_hat_net})
df_predictions.head()
df.drop('y', axis=1, inplace=False).columns
coefficients = pd.DataFrame({'Features':df.drop('y', axis=1, inplace=False).columns})
coefficients['model_lin']= model_linear.coef_
coefficients['model_ridge']= model_ridge.coef_
coefficients['model_lasso']= model_lasso.coef_
coefficients['model_net']= model_net.coef_
coefficients
评价各模型RMSE
MSE_test = np.mean(np.square(df_predictions['y_test'] - df_predictions['y_hat_linear']))
RMSE_test = np.sqrt(MSE_test)
np.round(RMSE_test,3)
0.159
MSE_test = np.mean(np.square(df_predictions['y_test'] - df_predictions['y_hat_ridge']))
RMSE_test = np.sqrt(MSE_test)
np.round(RMSE_test,3)
0.094
MSE_test = np.mean(np.square(df_predictions['y_test'] - df_predictions['y_hat_net']))
RMSE_test = np.sqrt(MSE_test)
np.round(RMSE_test,3)
1.048
绘制回归系数与 alpha 的对比图
1) 岭回归系数与阿尔法
alpha_ridge = 10**np.linspace(-2,4,100)
plt.figure(figsize=(12,5))
plt.plot(alpha_ridge,'or' )
plt.xlabel('steps')
plt.ylabel('alpha (lambda)')
plt.show()
ridge = Ridge()
coefs_ridge = []
for i in alpha_ridge:
ridge.set_params(alpha = i)
ridge.fit(X_train, y_train)
coefs_ridge.append(ridge.coef_)
np.shape(coefs_ridge)
(100, 13)
plt.figure(figsize=(12,10))
ax = plt.gca()
ax.plot(alpha_ridge, coefs_ridge)
ax.set_xscale('log')
plt.axis('tight')
plt.xlabel('alpha')
plt.ylabel('weights: scaled coefficients')
plt.title('Ridge regression coefficients Vs. alpha')
plt.legend(df.drop('y',axis=1, inplace=False).columns)
plt.show()
2) 套索回归系数与阿尔法
alpha_lasso = 10**np.linspace(-3,1,100)
lasso = Lasso()
coefs_lasso = []
for i in alpha_lasso:
lasso.set_params(alpha = i)
lasso.fit(X_train, y_train)
coefs_lasso.append(lasso.coef_)
np.shape(coefs_lasso)
(100, 13)
plt.figure(figsize=(12,10))
ax = plt.gca()
ax.plot(alpha_lasso, coefs_lasso)
ax.set_xscale('log')
plt.axis('tight')
plt.xlabel('alpha')
plt.ylabel('weights: scaled coefficients')
plt.title('Lasso regression coefficients Vs. alpha')
plt.legend(df.drop('y',axis=1, inplace=False).columns)
plt.show()
3) 弹性净回归系数与 alph
alpha_elasticnet = 10**np.linspace(-3,2,100)
elasticnet = ElasticNet()
coefs_elasticnet = []
for i in alpha_elasticnet:
elasticnet.set_params(alpha = i)
elasticnet.fit(X_train, y_train)
coefs_elasticnet.append(elasticnet.coef_)
np.shape(coefs_elasticnet)
(100, 13)
plt.figure(figsize=(12,10))
ax = plt.gca()
ax.plot(alpha_elasticnet, coefs_elasticnet)
ax.set_xscale('log')
plt.axis('tight')
plt.xlabel('alpha')
plt.ylabel('weights: scaled coefficients')
plt.title('Elastic Net regression coefficients Vs. alpha')
plt.legend(df.drop('y',axis=1, inplace=False).columns)
plt.show()
交叉验证(找lambda)
1) 岭回归的最优阿尔法
ridgecv = RidgeCV()
ridgecv.fit(X_train, y_train)
ridgecv.alpha_
0.1
alpha_ridge_opt = ridgecv.alpha_
2) 套索回归的最优阿尔法
lassocv = LassoCV()
lassocv.fit(X_train, y_train)
lassocv.alpha_
alpha_lasso_opt = lassocv.alpha_
0.0013485478432168933
3) 弹性净回归的最优 alpha
elasticnetcv = ElasticNetCV()
elasticnetcv.fit(X_train, y_train)
elasticnetcv.alpha_
0.0019027476880882372
elasticnetcv.l1_ratio
0.5
alpha_elasticnet_opt = elasticnetcv.alpha_
#可视化
from sklearn.model_selection import cross_val_score
import sklearn.metrics
RMSE_CV=[]
iterator=np.arange(0.0,0.02,0.001)
for i in iterator:
MSE=-cross_val_score(estimator=ElasticNet(alpha=i),X=X_train,y=y_train,cv=5,scoring="neg_mean_squared_error")
RMSE_CV.append(np.sqrt(MSE).mean())
output=pd.DataFrame(list(iterator),columns=['lambda_ElasticNet'])
output['RMSE_CV']=RMSE_CV
output.head()
output['RMSE_CV'].idxmin()
3
output['lambda_ElasticNet'][output['RMSE_CV'].idxmin()]
0.003
sns.lineplot(x='lambda_ElasticNet', y='RMSE_CV', data=output , color='r', label="RMSE_CV vs lambda_ridge")
plt.show()