Machine Learning #1 Linear Regression : 근로자 임금 회귀분석

 

GitHub - SeungbeomDo/DataAnalysis: Practical Codes for Data Analysis using Machine Learning and Deep Learning

---

1. 선형회귀모형의 도입

• 머신러닝 모델이라고 하면 복잡한 방법론을 떠올린다. 하지만 학부 통계학 강의에서도 쉽게 다루는 선형회귀모델도 머신러닝 모델의 한 유형이다. 모델이 주어진 데이터를 학습해서 예측의 구조를 정비하고 성능을 향상한다는 머신러닝의 일반적인 정의를 충족하기 때문이다.

 

• 선형회귀란 종속변수 $Y$와 $X$ 간의 선형적인 관계를 구하는 것을 말한다. 종속변수인 $Y$를 연구자가 선택한 k개의 독립변수 $X_{1}, X_{2}, ... , X_{k}$로 설명하는 선형회귀모형을 생각하자. 모형은 다음의 방정식 형태를 갖는다.

$$Y_{i} = \beta_{0} + \beta_{1}X_{i,1} + \beta_{2}X_{i,2} + ... + \beta_{k}X_{i,k} + \epsilon_{i}$$

 

• $Y_{i}$는 샘플 중 i번째 종속변수 관측값을 말한다. $\epsilon_{i}$는 모형으로 설명할 수 없는 잔차를 말한다. 
• 우리가 직접 추정해야 하는 대상은 독립변수들의 계수 집합 $\{\beta_{0}, \beta_{1}, ... , \beta_{k}\}$이다. 

 

• 2022년 경제활동인구조사 자료를 통해 근로자의 임금을 설명하는 선형회귀모형을 만들어보자. 임금과 상관이 높을 것으로 생각되는 성별, 학력, 기업규모, 정규직여부, 경력, 나이, 노조보호여부를 독립변수로 사용한다.

$$임금_{i} = \beta_{0} +\beta_{1}성별_{i} +\beta_{2}학력_{i} + \beta_{3}기업규모_{i} +\beta_{4}정규직여부_{i} +\beta_{5}경력_{i} +\beta_{6}나이_{i} +\beta_{7}노조보호여부_{i} + \epsilon_{i}$$

---

2. 데이터 전처리

2.1. 라이브러리 및 데이터 임포트

• 경제활동인구조사 데이터를 불러온 후, 근로자만 추리기 위해 임금이 0이 아닌 샘플들만 추출한다.

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns

data = pd.read_csv('path')
data = data[data['income'] > 0]

 

2.2. 변수 처리

• 성별: 남성은 1, 여성은 2의 값을 갖는 범주형 변수로 라벨링되어 있으므로, 남성 여부를 더미화한다.
• 정규직 여부: 상용직은 1, 임시직은 2, 일용직은 3으로 범주형 라벨링 돼있으므로, 정규직 여부로 더미화한다.
• 노조 보호 여부: 노조가 없는 경우는 1, 노조 가입한 경우는 2, 노조 가입 대상이지만 가입하지 않은 경우는 3, 노조가 있지만 가입대상이 아닌 경우는 4로 라벨링 돼있으므로 노조 보호 대상(2, 3) 여부로 더미화한다.
• 경력: 2022년 - 근로 시작 연도로 처리한다.

#성별
data['gender'] = (data['gender']==1)

#정규직 여부
data['regular'] = (data['worktype'] == 1)

#노조 보호 여부
data['union'] = (data['union'] == 2) | (data['union'] == 3)

#경력
data['workstart'] = data['workstart'].astype('str')
data['workstart'] = [x[:4] for x in data['workstart']]
data['workstart'] = data['workstart'].astype('int')
data['exp'] = 2022 - data['workstart']
data = data.drop(columns = 'workstart')

---

3. 기술분석을 통한 독립변수 점검

• 사용할 독립변수를 다시 한 번 점검할 필요가 있다. 선형회귀모형에서 (1) 설명력이 없는 독립변수를 포함하거나 (2) 다른 독립변수들에 의해 설명되는 독립변수를 사용하는 경우 일정한 문제가 발생한다.

 

• (1) 설명력이 없는 독립변수는 모델에 noise로 작용하여 학습 효과를 저해시킨다. 또한 독립변수의 개수가 많을수록 과대적합의 위험이 높아지고 계산비용이 커지는 차원의 저주 문제가 나타나는데, 이를 막기 위해서는 불필요한 독립변수는 걸러내어 모형의 차원을 줄여주어야 한다.
• (2) 다른 독립변수들에 의해 설명되는 독립변수가 존재할 경우, 다중공선성(multicollinearity)가 존재한다고 한다. 다중공선성의 가장 큰 문제는 추정된 회귀계수의 해석을 부정확하게 만든다는 것이다. 회귀계수는 다른 변수들이 고정될 때 주어진 변수의 변화가 종속변수에 미치는 영향력을 측정한다. 그런데 복수의 독립변수들이 상관된 움직임을 보인다면 이런 식의 해석이 불가능하다. 또한 다중공선성은 추정된 회귀계수의 표준오차를 크게 만들어 통계적 유의성을 검정하기 어렵게 만든다는 문제도 있다.

 

• 언급한 문제들을 예방하기 위해서는 1차적으로 기술분석을 통해 변수 간 상관계수를 파악하는 것이다. (1) 종속변수와 상관이 너무 낮은 독립변수는 불필요한 독립변수라고 볼 수 있고 (2) 다른 독립변수와 상관이 너무 큰 독립변수는 다중공선성을 유발하는 원인이 된다.

data.corr()

 

• 종속변수인 임금과 특별히 상관계수가 너무 낮거나, 혹은 다른 변수들과 상관관계가 유의하게 높은 독립변수는 보이지 않는다. 따라서 별다른 처리 없이 진행하도록 한다.

---

4. 훈련 데이터셋 분리

• train set과 test set을 분리하여 정확한 모형의 성능을 측정하도록 한다.

X = data.iloc[:,:-1] #income 제외
y = data.iloc[:,-1] #income만 포함

from sklearn.model_selection import train_test_split

train_X, test_X, train_y, test_y = train_test_split(X, y, train_size=0.8,
                                                    test_size=0.2, random_state=1)
print(train_X.shape, test_X.shape, train_y.shape, test_y.shape)

---

5. 과대적합 예방: 규제항 포함

• 선형회귀모형의 과대적합을 방지하는 주요한 방법은 독립변수 벡터의 차원을 줄여주는 것이다. 이때 독립변수 벡터의 차원을 줄이기 위한 테크닉으로 규제항을 사용할 수 있다.
• 일반적으로 선형회귀모델을 훈련할 때는 예측오차의 제곱합이 비용함수가 된다. 규제항을 포함하는 것은, 모델 훈련 시 개별 독립변수들의 회귀계수 합이 너무 커지지 않도록 비용함수를 수정하는 것을 의미한다. 예측오차를 줄이면서 동시에 회귀계수 합도 줄이기 위해 불필요한 독립변수들의 회귀계수는 0에 가까워지도록 훈련하게 된다.

 

•   대표적인 규제 방법은 Ridge 회귀와 Lasso 회귀가 있다. 두 규제 회귀의 비용함수는 아래와 같다.

 

Ridge: $\Sigma_{i=1}^{n} (y_{i} - \hat{y_{i}})^{2} + \lambda \Sigma_{j=1}^{k}\hat{\beta}_{j}^{2}$

Lasso: $\Sigma_{i=1}^{n} (y_{i} - \hat{y_{i}})^{2} + \lambda \Sigma_{j=1}^{k}|\hat{\beta}_{j}|$

 

• 두 규제 방법의 차이는 Ridge 회귀의 경우 불필요한 변수의 회귀계수를 0으로 만들 수는 없지만 Lasso의 경우는 그것도 가능하다는 데 있다.
• 한편 규제항에 곱해진 파라미터 $\lambda$는 규제의 강도를 측정한다. 클수록 규제가 강해진다.

---

6. 모델 훈련 및 검증

6.1. 훈련 방법

• 선형회귀모델을 훈련한다는 것은 train set으로부터 예측오차 제곱합을 최소화(OLS)하는 베타값을 찾는다는 것을 말한다. 예측오차 제곱합을 최소화하는 방법은 (1) 1계 조건(First order Condition)을 충족하는 해석적 해를 찾거나, (2) 해석적 해가 존재하지 않는 경우 유사역행렬을 사용하거나, (3) 비용함수를 최소화하는 파라미터 값을 점진적으로 찾아가는 수치적 방법(경사하강법)으로 세분화할 수 있다.

 

6.2. 파라미터 튜닝

• 한편 규제항 회귀를 사용할 것인데, 파라미터 $\lambda$의 최적값이 얼마인지 모르니 후보값들을 두고 $R^{2}$값을 최대화하는 값을 찾아보도록 한다. 이런 과정을 파라미터 튜닝이라고도 한다.

 

6.3. 검증 데이터셋(vallidation set)의 활용

• 훈련과정에서 test set을 건드리지 않고 좋은 파라미터를 찾기 위해, 그리고 트레인셋만을 활용하면서도 트레인셋에 과대적합되는 것을 막기 위해 검증 데이터셋을 활용한다.
• train set에서 다시 subset으로 검증 데이터셋을 만들어서 이를 제외한 데이터로 훈련을 하고 검증 데이터셋으로 성능을 측정하는 것을 말한다. 
• 주로 k-fold 방법을 사용하는데 훈련 데이터셋에서 임의로 검증 데이터셋을 뽑아서 훈련을 해보고, 그 다음에 다시 임의로 새로운 검증 데이터셋을 뽑아서 훈련을 해보고, ... , 를 k번 반복한다는 뜻이다.

 

• 가장 좋은 규제 파라미터를 결정하기 위해 k-fold 방법을 적용해보았다.

from sklearn.model_selection import cross_val_score
from sklearn.linear_model import Ridge, Lasso

alpha_list = [0.001, 0.01, 0.1, 1, 10]

from sklearn import metrics
from sklearn.metrics import r2_score # R^2

ridge_r2_list_train = []
ridge_r2_list_test = []
lasso_r2_list_train = []
lasso_r2_list_test = []

import warnings
warnings.filterwarnings('ignore')

for alpha in alpha_list:

    ridge_model = Ridge(alpha=alpha)
    lasso_model = Lasso(alpha=alpha)

    ridge_r2_list_train.append(np.mean(cross_val_score(ridge_model,train_X,train_y,scoring='r2',cv=10)))
    lasso_r2_list_train.append(np.mean(cross_val_score(lasso_model,train_X,train_y,scoring='r2',cv=10)))

plt.plot(alpha_list, ridge_r2_list_train, label='Ridge R2 Train')
plt.plot(alpha_list, lasso_r2_list_train, label='Lasso R2 Train')

plt.legend()

• 규제 파라미터의 후보를 0.001, 0.01, 0.1, 1, 10으로 두고 Ridge 회귀와 Lasso 회귀를 모두 실험해보았는데 둘 다 파라미터가 0.001일 때 10개 검증 데이터셋의 $R^{2}$ 평균이 가장 높게 나왔다. 왜냐면 아까 상관계수 행렬을 통해 고찰하였듯이, 독립변수들이 다 유의미해서 굳이 규제항을 포함할 필요가 없기 때문이다.

 

• 따라서 그냥 규제 회귀를 실시하지 않고 일반적인 선형회귀를 돌렸다. sklearn 라이브러리의 선형회귀모델을 사용할 경우 경사하강법으로 회귀계수의 수치적 해를 구해준다고 한다.

from sklearn.metrics import mean_squared_error, mean_absolute_error

def evaluate_regr(y,pred):
    rmse_val = np.sqrt(mean_squared_error(y,pred))
    mae_val = mean_absolute_error(y,pred)

    print('RMSE: {0:.3F}, MAE: {1:.3F}'.format(rmse_val, mae_val))
    return (rmse_val, mae_val)
    
    
from sklearn.linear_model import LinearRegression

lr_model = LinearRegression()
lr_model.fit(train_X, train_y) 
pred = lr_model.predict(test_X) 

result = evaluate_regr(test_y, pred)
plt.bar(['RMSE','MAE'], result)

---

7. 모형 해석

• 전통적인 선형회귀모델이 우수한 이유 중 하나는, 모형을 설명하는 것이 가능하다는 점이다. 모형에서 추정된 회귀계수들을 알고 있기 때문에 어떤 변수가 중요한지 알아낼 수 있다. 계수값이 크면 중요하다.

coef = lr_model.coef_
coef_abs = np.abs(coef)

idx = coef_abs.argsort()[::-1]

coef_abs = coef_abs[idx]

name = X.head(1).columns[idx].values

sns.barplot(name[0:10], coef_abs[0:10])
plt.xticks(rotation=45)
plt.ylabel('Reg_Coef')
plt.title('Reg_Coef per features; TOP10')
plt.show()

 

• 참 애석하게도, 성별과 정규직 여부따위가 임금을 결정할 때 제일 중요하다고 나온다. 물론 여러 설명변수들을 다 포함한 것도 아니고, 이론적 고찰 없이 내 맘대로 한 거라 이런 결론이 절대적인 것은 결코 아니다.

---

8. 스케일링

8.1. 스케일링의 의의

• 스케일링(Scaling)이란 변수들의 단위들을 표준화해주는 것을 말한다. 

 

• 예측과제를 수행할 때는 데이터의 원 단위를 유지하는 것이 좋겠지만 변수 중요도를 해석하는 것이 메인 과제라면 스케일링이 필요하다. 각 변수마다 단위 차이가 크다면 이 점이 계수의 추정에도 영향을 미치기 때문이다.
  
  • 예시를 들자. 어떤 사람의 몸무게(kg)를 설명하는 데 독립변수로 일일 평균 섭취 칼로리(kcal)와 키(cm)를 사용한다고 하자. 칼로리는 기본적으로 1천 단위의 값을 갖지만 키는 기껏해야 180, 190 정도이다. 이때 100이 조금 안 되는 단위를 가진 몸무게를 맞추려면 칼로리에 붙는 베타값은 매우 작아야 할 것이다. 반면 키에는 상대적으로 큰 베타값이 붙게 될 것이다. 
  • 이런 경우 몸무게를 예측하는 것에는 문제가 안 된다. 하지만 변수 중요도 해석의 맥락에서 보자. 칼로리에 붙는 베타값이 키에 붙는 베타값보다 작다고 해서, 칼로리가 덜 중요하다고 볼 수 있을까? 이런 문제를 해결하기 위해 스케일링을 사용한다.

 

8.2. Min-Max Scaling

• 스케일링에는 대표적으로 두 가지가 있는데 (1) min-max scaling과 (2) standard scaling이다.

 

• (1) min-max scaling: j번째 독립변수의 i번째 관측값 $X_{i,j}$에 대하여 min-max 스케일링은 다음과 같다.

$$\frac{max(X_{j})-X_{i,j}}{max(X_{j})-min(X_{j})}$$

 

• 이때 스케일된 값은 0과 1 사이의 범위를 갖는다. 당연히 더미변수는 빼놓고 해야 한다.

from sklearn.preprocessing import MinMaxScaler

scale_data = data.drop(columns = ['gender', 'regular', 'union'])
remaining_data = data[['gender', 'regular', 'union']]

scaler = MinMaxScaler()
scaler.fit(scale_data)
trans_data = scaler.transform(scale_data)

trans_data_df = pd.DataFrame(data=trans_data, columns=scale_data.columns)
trans_data_df = pd.concat([trans_data_df.reset_index(drop=True), remaining_data.reset_index(drop=True)], axis = 1)
trans_data_df

 

• 베타값의 해석은 해당 변수가 1%p 증가했을 때, 종속변수에 미치는 %p 영향력이다.

• 원 단위로 해석했을 때와는 반대로, 성별과 정규직 여부의 중요성은 감소하고 경력과 학력 등 합리적 요인의 중요도가 높아졌다. 세상 참 다행이다...

 

8.3. Standard Scaling

• (2) min-max scaling: j번째 독립변수의 i번째 관측값 $X_{i,j}$에 대하여 standard 스케일링은 다음과 같다.

$$\frac{X_{i,j}-E[X_{j}]}{std(X_{j})}$$

 

• 이때 스케일된 값은 평균이 0이고 표준편차가 1인 임의의 분포를 따르게 된다.
• 베타값의 해석은 해당 변수가 1 표준편차 증가했을 때, 종속변수에 미치는 표준편차 영향력이다.

• 이 경우에는 다시 성별과 정규직 여부의 중요도가 높아졌다.