kaggle 주택 가격 예측(3) - 간단한 regression으로 예측하기(상위 30%)

 

저번 글에는 데이터를 전처리하고 저장했었다면 이번 글에서는 전처리 데이터로 학습하고 kaggle에 제출해 볼 예정이다. 복잡한 알고리즘 이전에 간단한 머신러닝으로 돌린 성능을 확인하기 위해, 이번 글에서는 우선 Linear regression, ridge regression, rasso regression, elasticnet 네 가지로 실습을 해보자.

 

그나저나 노트북을 그대로 옮겨오면 사이즈 때문에 별로 예쁘지도 않고, 작성 후 재수정시, html구조가 깨져 노트북 테마가 다 사라져 버린다. 다른 방법을 강구해봐야겠다.

 

---

 

 

 

 

 

4 kinds of regression¶

In [1]:

# Imports
import pandas as pd
import numpy as np
from sklearn.model_selection import cross_val_score, train_test_split, KFold
from sklearn.preprocessing import StandardScaler, RobustScaler
from sklearn.linear_model import LinearRegression, RidgeCV, LassoCV, ElasticNetCV
from sklearn.metrics import mean_squared_error, make_scorer
from IPython.display import display
import matplotlib.pyplot as plt

import warnings
def ignore_warn(*args, **kwargs): pass
warnings.warn = ignore_warn #ignore annoying warning (from sklearn and seaborn)

%matplotlib inline

pd.set_option('display.float_format', lambda x: '%.3f' % x)
pd.set_option('display.max_colwidth', -1)
pd.set_option('display.max_rows', 500)

 

In [2]:

# 데이터 읽어오기
# 이전 글에서 데이터 전처리를 잘 따라왔다면 해당경로에 데이터가 저장되어 있을것이다.
# 동일한 데이터가 아니더라도 똑같은 포맷이라면 동일한 방식으로 사용해 볼 수 있다.
X_train = pd.read_csv("../data/X_train.csv", index_col=0, header=0)
X_test = pd.read_csv("../data/X_test.csv", index_col=0, header=0)
y_train = pd.read_csv("../data/y_train.csv", header=None, index_col = 0).values[:, 0]
Id_test = pd.read_csv("../data/Id_train.csv",  index_col=0, header=None).values[:, 0]
print(X_train.shape)
print(X_test.shape)
print(y_train.shape)
print(Id_test.shape)

 

(1458, 333)
(1459, 333)
(1458,)
(1459,)

 

Modeling

 

In [3]:

# 스케일러를 사용해주어야 오버플로우나 언더플로우를 방지할 수있다.
# 스케일러를 사용 하기전에는 아웃라이어가 무조건 제거되어야한다.
# 그런데 나의 경우 아웃라이어를 심도깊게 보지 않았으므로 그나마 아웃라이어에 강한 Robust Scaler를 사용한다.
# 모델내에서 normalize = True를 주는 방법도 있다.(이게 더 잘나오는 경우도 더 많다.)
categorical_features = X_train.select_dtypes(include = ["object"]).columns
numerical_features = X_train.select_dtypes(exclude = ["object"]).columns

stdSc = StandardScaler().fit(X_train.loc[:, numerical_features])
X_train.loc[:, numerical_features] = stdSc.transform(X_train.loc[:, numerical_features])
X_test.loc[:, numerical_features] = stdSc.transform(X_test.loc[:, numerical_features])

 

In [4]:

n_folds = 10
def rmse_cv(model):
    kf = KFold(n_folds, shuffle=True, random_state=42).get_n_splits(X_train.values)
    rmse= np.sqrt(-cross_val_score(model, X_train.values, y_train, scoring="neg_mean_squared_error", cv = kf))
    return(rmse)

 

1* Linear Regression without regularization

 

In [5]:

# Linear Regression
lr = LinearRegression()
lr.fit(X_train, y_train)

# Look at predictions on training and validation set
print("RMSE on Training set :", rmse_cv(lr).mean())
y_pred = lr.predict(X_train)

# Plot residuals
plt.scatter(y_pred, y_pred - y_train, c = "blue", marker = "s", label = "Training data")
plt.title("Linear regression")
plt.xlabel("Predicted values")
plt.ylabel("Residuals")
plt.legend(loc = "upper left")
plt.hlines(y = 0, xmin = 10.5, xmax = 13.5, color = "red")
plt.show()

# Plot predictions
plt.scatter(y_pred, y_train, c = "blue", marker = "s", label = "Training data")
plt.title("Linear regression")
plt.xlabel("Predicted values")
plt.ylabel("Real values")
plt.legend(loc = "upper left")
plt.plot([10.5, 13.5], [10.5, 13.5], c = "red")
plt.show()

 

RMSE on Training set : 126649931115.9172

 

 

 

RMSE가 꽤나 높다. Linear Regression은 이 데이터에 맞지 않다고 판단, 사용하지 않는것으로 결론 내리자.

 

2* Linear Regression with Ridge regularization (L2 penalty)

Ridge regression 은 L2 패널티를 주는 리니어 모델로 cost function의 제곱에 weight를 곱하는 방식으로 작동한다.

 

In [6]:

# 2* Ridge
ridge = RidgeCV(alphas = [0.01, 0.03, 0.06, 0.1, 0.3, 0.6, 1, 3, 6, 10, 15, 25, 30, 35, 40, 45, 50, 55, 60])
ridge.fit(X_train, y_train)
alpha = ridge.alpha_
print("Best alpha :", alpha)

print("Try again for more precision with alphas centered around " + str(alpha))
ridge = RidgeCV(alphas = [alpha * .6, alpha * .65, alpha * .7, alpha * .75, alpha * .8, alpha * .85, 
                          alpha * .9, alpha * .95, alpha, alpha * 1.05, alpha * 1.1, alpha * 1.15,
                          alpha * 1.25, alpha * 1.3, alpha * 1.35, alpha * 1.4], cv = 10)
ridge.fit(X_train, y_train)
alpha = ridge.alpha_
print("Best alpha :", alpha)

print("Ridge RMSE on Training set :", rmse_cv(ridge).mean())
y_rdg = ridge.predict(X_train)

# Plot residuals
plt.scatter(y_rdg, y_rdg - y_train, c = "blue", marker = "s", label = "Training data")
plt.title("Linear regression with Ridge regularization")
plt.xlabel("Predicted values")
plt.ylabel("Residuals")
plt.legend(loc = "upper left")
plt.hlines(y = 0, xmin = 10.5, xmax = 13.5, color = "red")
plt.show()

# Plot predictions
plt.scatter(y_rdg, y_train, c = "blue", marker = "s", label = "Training data")
plt.title("Linear regression with Ridge regularization")
plt.xlabel("Predicted values")
plt.ylabel("Real values")
plt.legend(loc = "upper left")
plt.plot([10.5, 13.5], [10.5, 13.5], c = "red")
plt.show()

# Plot important coefficients
coefs = pd.Series(ridge.coef_, index = X_train.columns)
print("Ridge picked " + str(sum(coefs != 0)) + " features and eliminated the other " +  \
      str(sum(coefs == 0)) + " features")
imp_coefs = pd.concat([coefs.sort_values().head(10),
                     coefs.sort_values().tail(10)])
imp_coefs.plot(kind = "barh")
plt.title("Coefficients in the Ridge Model")
plt.show()

 

Best alpha : 60.0
Try again for more precision with alphas centered around 60.0
Best alpha : 84.0
Ridge RMSE on Training set : 0.113967023647391

 

 

 

Ridge picked 332 features and eliminated the other 1 features

 

 

기본 Linear Regression 보다는 훨씬 낫다. 확실히 regularization이 오버핏을 잘 방지했나보다. 그런데 alpha 값이 조금 의하하다.

 

3* Linear Regression with Lasso regularization (L1 penalty)

LASSO 은 다른 정규화 방식을 사용하는데 일정수치 이하로 적은 cost를 무시하는 방향으로 작용한다.

 

In [7]:

# 3* Lasso
lasso = LassoCV(alphas = [0.0001, 0.0003, 0.0006, 0.001, 0.003, 0.006, 0.01, 0.03, 0.06, 0.1, 
                          0.3, 0.6, 1], 
                max_iter = 50000, cv = 10, n_jobs = -1)
lasso.fit(X_train, y_train)
alpha = lasso.alpha_
print("Best alpha :", alpha)

print("Try again for more precision with alphas centered around " + str(alpha))
lasso = LassoCV(alphas = [alpha * .6, alpha * .65, alpha * .7, alpha * .75, alpha * .8, 
                          alpha * .85, alpha * .9, alpha * .95, alpha, alpha * 1.05, 
                          alpha * 1.1, alpha * 1.15, alpha * 1.25, alpha * 1.3, alpha * 1.35, 
                          alpha * 1.4], 
                max_iter = 50000, cv = 10, n_jobs=-1)
lasso.fit(X_train, y_train)
alpha = lasso.alpha_
print("Best alpha :", alpha)

print("Lasso RMSE on Training set :", rmse_cv(lasso).mean())
y_las = lasso.predict(X_train)

# Plot residuals
plt.scatter(y_las, y_las - y_train, c = "blue", marker = "s", label = "Training data")
plt.title("Linear regression with Lasso regularization")
plt.xlabel("Predicted values")
plt.ylabel("Residuals")
plt.legend(loc = "upper left")
plt.hlines(y = 0, xmin = 10.5, xmax = 13.5, color = "red")
plt.show()

# Plot predictions
plt.scatter(y_las, y_train, c = "blue", marker = "s", label = "Training data")
plt.title("Linear regression with Lasso regularization")
plt.xlabel("Predicted values")
plt.ylabel("Real values")
plt.legend(loc = "upper left")
plt.plot([10.5, 13.5], [10.5, 13.5], c = "red")
plt.show()

# Plot important coefficients
coefs = pd.Series(lasso.coef_, index = X_train.columns)
print("Lasso picked " + str(sum(coefs != 0)) + " features and eliminated the other " +  \
      str(sum(coefs == 0)) + " features")
imp_coefs = pd.concat([coefs.sort_values().head(10),
                     coefs.sort_values().tail(10)])
imp_coefs.plot(kind = "barh")
plt.title("Coefficients in the Lasso Model")
plt.show()

 

Best alpha : 0.003
Try again for more precision with alphas centered around 0.003
Best alpha : 0.0033000000000000004
Lasso RMSE on Training set : 0.10761950086364672

 

 

 

Lasso picked 102 features and eliminated the other 231 features

 

 

RMSE 결과가 살짝 더 좋다. 흥미로운 점은 231개의 feature 를 버리고 100개 정도의 feature만 사용했다는 것이다. 나중에 필요 없었던 feature를 확인 해 보도록 하자.

 

4* Linear Regression with ElasticNet regularization (L1 and L2 penalty)

ElasticNet 은 Ridge와 Lasso의 절충점이다. L1 패널티와 L2 패널티 모두를 고려한다.

 

In [8]:

# 4* ElasticNet

elasticNet = ElasticNetCV(l1_ratio = [0.1, 0.3, 0.5, 0.6, 0.7, 0.8, 0.85, 0.9, 0.95, 1],
                          alphas = [0.0001, 0.0003, 0.0006, 0.001, 0.003, 0.006, 
                                    0.01, 0.03, 0.06, 0.1, 0.3, 0.6, 1, 3, 6], 
                          max_iter = 50000, cv = 10, n_jobs=-1)
elasticNet.fit(X_train, y_train)
alpha = elasticNet.alpha_
ratio = elasticNet.l1_ratio_
print("Best l1_ratio :", ratio)
print("Best alpha :", alpha )
print("Try again for more precision with l1_ratio centered around " + str(ratio))
elasticNet = ElasticNetCV(l1_ratio = [ratio * .85, ratio * .9, ratio * .95, ratio, ratio * 1.05, ratio * 1.1, ratio * 1.15],
                          alphas = [0.0001, 0.0003, 0.0006, 0.001, 0.003, 0.006, 0.01, 0.03, 0.06, 0.1, 0.3, 0.6, 1, 3, 6], 
                          max_iter = 50000, cv = 10, n_jobs=-1)
elasticNet.fit(X_train, y_train)
if (elasticNet.l1_ratio_ > 1):
    elasticNet.l1_ratio_ = 1    
alpha = elasticNet.alpha_
ratio = elasticNet.l1_ratio_
print("Best l1_ratio :", ratio)
print("Best alpha :", alpha )

print("Now try again for more precision on alpha, with l1_ratio fixed at " + str(ratio) + 
      " and alpha centered around " + str(alpha))
elasticNet = ElasticNetCV(l1_ratio = ratio,
                          alphas = [alpha * .6, alpha * .65, alpha * .7, alpha * .75, alpha * .8, alpha * .85, alpha * .9, 
                                    alpha * .95, alpha, alpha * 1.05, alpha * 1.1, alpha * 1.15, alpha * 1.25, alpha * 1.3, 
                                    alpha * 1.35, alpha * 1.4], 
                          max_iter = 50000, cv = 10, n_jobs=-1)
elasticNet.fit(X_train, y_train)
if (elasticNet.l1_ratio_ > 1):
    elasticNet.l1_ratio_ = 1    
alpha = elasticNet.alpha_
ratio = elasticNet.l1_ratio_
print("Best l1_ratio :", ratio)
print("Best alpha :", alpha )

print("ElasticNet RMSE on Training set :", rmse_cv(elasticNet).mean())
y_ela = elasticNet.predict(X_train)

# Plot residuals
plt.scatter(y_ela, y_ela - y_train, c = "blue", marker = "s", label = "Training data")
plt.title("Linear regression with ElasticNet regularization")
plt.xlabel("Predicted values")
plt.ylabel("Residuals")
plt.legend(loc = "upper left")
plt.hlines(y = 0, xmin = 10.5, xmax = 13.5, color = "red")
plt.show()

# Plot predictions
plt.scatter(y_train, y_ela, c = "blue", marker = "s", label = "Training data")
plt.title("Linear regression with ElasticNet regularization")
plt.xlabel("Predicted values")
plt.ylabel("Real values")
plt.legend(loc = "upper left")
plt.plot([10.5, 13.5], [10.5, 13.5], c = "red")
plt.show()

# Plot important coefficients
coefs = pd.Series(elasticNet.coef_, index = X_train.columns)
print("ElasticNet picked " + str(sum(coefs != 0)) + " features and eliminated the other " +  str(sum(coefs == 0)) + " features")
imp_coefs = pd.concat([coefs.sort_values().head(10),
                     coefs.sort_values().tail(10)])
imp_coefs.plot(kind = "barh")
plt.title("Coefficients in the ElasticNet Model")
plt.show()

 

Best l1_ratio : 1.0
Best alpha : 0.003
Try again for more precision with l1_ratio centered around 1.0
Best l1_ratio : 1
Best alpha : 0.003
Now try again for more precision on alpha, with l1_ratio fixed at 1 and alpha centered around 0.003
Best l1_ratio : 1
Best alpha : 0.0033000000000000004
ElasticNet RMSE on Training set : 0.10761950086364672

 

 

 

ElasticNet picked 102 features and eliminated the other 231 features

 

 

ElasticNet은 L1 패널티만 주었을때와 거의 비슷한 결과를 나타낸다. 즉 lasso regression 과 성능이 같은 샘이다.

 

In [9]:

pred = elasticNet.predict(X_test)
pd_test_pred = pd.DataFrame({'Id': Id_test, 'SalePrice': np.expm1(pred)})

 

테스트 데이터로도 예측한 결과를 제출해보자.
rmse 0.13204으로 상위 35%정도에 랭크된다. 만족스럽진 못한 결과다.

 

In [10]:

pd_test_pred.to_csv('submission.csv', index=False)

 

결론

캐글에서 적당한 노트북을 발견해서 regression 의 대표적인 4가지 알고리즘을 사용해보았다.
간단한 elasticNet으로 제출시 rmse 0.13204으로 상위 35% 정도에 랭크된다.
다음 글에는 xgboost 같은 더 유명한 알고리즘을 사용하고, Stack 기법과 ensemble 기법을 이용해 더 나은 결과물을 만들어 내보자.

 

In [11]:

# 티스토리에 올리기 위한 레이아웃 조정
from IPython.core.display import display, HTML
display(HTML("<style>.container {width:95% !important;}</style>"))

 

 

---

 

l1, l2 penalty 모두를 고려한 ElasticNet으로 제출햇을때 상위 35% 정도에 랭크된다(결과는 사실상 L2 페널티는 없는 것과 같다).

썩 만족스러운 결과는 아닌데 다음 글에서는 요즘 캐글에서 핫한 알고리즘들을 써서 제출해 봐야겠다.