2021-09-18(따릉이 프로젝트 완성하기 7)

*주요 요약

1. hour_bef_ozone의 결측치 예측

2. 결측치 다수 존재하는 변수들은 미 사용하기로 판단함

 

 

* 본 내용

결측치가 많은 ozone, pm10, pm2.5를 예측하여 값을 채우려 시도하고 있습니다. 기존에는 각 시간, 변수 별 평균값으로 결측치를 채웠습니다. 하지만 이러한 방법은 정확성을 떨어뜨릴수 있습니다 게다가 그러한 결측치가 워낙 많다 보니까 새로운 방법을 찾아야 했습니다. ozone, pm10, pm2.5를 제외한 feature을 활용했습니다.(id, count 제외)

 

먼저, 사용될 X feature들의 결측치 먼저 처리해야 합니다. 결측치 수가 많지 않으므로 이는 기존에 만들어 사용했던 시간,변수별 평균값 함수를 사용하여 채웠습니다. 아래 step 1을 통해 상세하게 설명드립니다.

 

Step 1. nan_mean_fill 함수를 사용하여 X feature 채우기 (단, precipitation은 분류 이므로 이 또한 제외)

def nan_mean_fill(dataset):
  #결측치가 있는 컬럼
  nan_columns = [col for col in dataset.iloc[:,2:-1].columns.to_list() if col!='hour_bef_precipitation' and col!='hour_bef_ozone' and col!='hour_bef_pm10' and col!='hour_bef_pm2.5']

  #for 문을 통해 결측치를 각 시간대별 평균값으로 대체
  for col in nan_columns:
    hour_mean_value = dataset.groupby('hour').mean()[col]
    hour_mean_value.fillna(hour_mean_value.mean(),inplace=True) #회귀를 통해서 값 채우기 시도

    for nan_hour in dataset[dataset[col].isna()]['hour'].unique():
      
      #nan index 구하기
      index= dataset[dataset[col].isna()].loc[dataset['hour']==nan_hour,col].index

      #채우기 
      for idx in index:
        dataset.loc[idx, col] = hour_mean_value[nan_hour]
        
  print(dataset.isna().sum())
  return dataset

train = nan_mean_fill(train)
test= nan_mean_fill(test)

 

id                          0
hour                        0
hour_bef_temperature        0
hour_bef_precipitation      2
hour_bef_windspeed          0
hour_bef_humidity           0
hour_bef_visibility         0
hour_bef_ozone             76
hour_bef_pm10              90
hour_bef_pm2.5            117
count                       0
dtype: int64
id                         0
hour                       0
hour_bef_temperature       0
hour_bef_precipitation     1
hour_bef_windspeed         0
hour_bef_humidity          0
hour_bef_visibility        0
hour_bef_ozone            35
hour_bef_pm10             37
hour_bef_pm2.5            36

precipitation과 target으로 사용될 ozone,pm10, pm2.5를 제외한 나머지 변수들의 값을 채운걸 확인해 볼수 있습니다.

 

다음은 Logistic을 통해 precipitation의 결측치를 채우겠습니다. 참고로 precipitation은 비가 온다 안온다에 대한 데이터이며 0과1로 구성되어 있습니다.

 

 

 

step 2. precipitation 결측값 채우기

from sklearn.linear_model import LogisticRegression
from sklearn import metrics

#x,y 분리
precipitation_x_del = pd.concat([train,test], axis=0).iloc[:,1:-4].dropna()
cls_x_train = precipitation_x_del[[col for col in precipitation_x_del.columns if col!='hour_bef_precipitation']]
cls_y_train =precipitation_x_del['hour_bef_precipitation']

#data split
cls_X_train, cls_X_val, cls_Y_train, cls_Y_test = train_test_split(cls_x_train, cls_y_train, test_size = 0.3, shuffle = True, random_state = 2021)

# model fit
cls_model = LogisticRegression()
cls_model.fit(cls_X_train, cls_Y_train)

#predict
y_pred = pd.Series(cls_model.predict(cls_X_val))

# model evaluation
print("Accuracy:", metrics.accuracy_score(cls_Y_test, y_pred))

Accuracy: 0.9662576687116564

 

# 불필요한 변수 제거
hour_bef_precipitation_train_df = train[train['hour_bef_precipitation'].isna()].iloc[:,1:-4]
del hour_bef_precipitation_train_df['hour_bef_precipitation']
display(hour_bef_precipitation_train_df)

hour_bef_precipitation_test_df = test[test['hour_bef_precipitation'].isna()].iloc[:,1:-3]
del hour_bef_precipitation_test_df['hour_bef_precipitation']
display(hour_bef_precipitation_test_df)

#predict -> 비가 오지 않음(0) 
print(cls_model.predict(hour_bef_precipitation_train_df))
print(cls_model.predict(hour_bef_precipitation_test_df))

# 결측값 0으로 채우기
train['hour_bef_precipitation'].fillna({934:0, 1035:0}, inplace=True)
print(train.isna().sum())

test['hour_bef_precipitation'].fillna({653:0}, inplace=True)
print(test.isna().sum())

#predict -> 비가 오지 않음(0) 
print('train의 hour_bef_precipitation 예측값:',cls_model.predict(hour_bef_precipitation_train_df))
print("test의 hour_bef_precipitation 예측값:",cls_model.predict(hour_bef_precipitation_test_df))

 

	hour	hour_bef_temperature	hour_bef_windspeed	hour_bef_humidity	hour_bef_visibility
934	0	14.788136	1.965517	58.169492	1434.220339
1035	18	20.926667	3.838333	40.450000	1581.850000
hour	hour_bef_temperature	hour_bef_windspeed	hour_bef_humidity	hour_bef_visibility
653	19	26.110345	3.541379	47.689655	1561.758621
[0. 0.]
[0.]
id                          0
hour                        0
hour_bef_temperature        0
hour_bef_precipitation      0
hour_bef_windspeed          0
hour_bef_humidity           0
hour_bef_visibility         0
hour_bef_ozone             76
hour_bef_pm10              90
hour_bef_pm2.5            117
count                       0
dtype: int64
id                         0
hour                       0
hour_bef_temperature       0
hour_bef_precipitation     0
hour_bef_windspeed         0
hour_bef_humidity          0
hour_bef_visibility        0
hour_bef_ozone            35
hour_bef_pm10             37
hour_bef_pm2.5            36
dtype: int64
train의 hour_bef_precipitation 예측값: [0. 0.]
test의 hour_bef_precipitation 예측값: [0.]

결측값을 채우기 위한 분류 결과, 모두 비가 오지 않았다는 예측을 얻었을 수가 있었습니다. 따라서 precipitation 결측값에 모두 0을 채웠습니다. 

 

X feature들의 결측값을 모두 채웠으므로 ozone, pm10, pm2.5 각각을 Y_target으로 설정하여 예측을 해보록 하겠습니다. 

 

step3. ozone, pm10, pm2.5 결측값 예측

def evaluate(y_true, y_pred):
  y_true, y_pred = np.array(y_true), np.array(y_pred)
  mape= np.mean(np.abs((y_true-y_pred)/y_true))*100
  mae= mean_absolute_error(y_true,y_pred)
  mse= mean_squared_error(y_true,y_pred)
  score = pd.DataFrame([mape, mae, mse], index=['mape','mae','mse'], columns=['score']).transpose()
  return score

 

def Nan_prediction(targets):
  # train 데이터 만들기
  for target in targets:
    if target =='hour_bef_pm2.5':
       use_columns = ['hour', 'hour_bef_temperature', 'hour_bef_precipitation','hour_bef_windspeed', 'hour_bef_humidity', 'hour_bef_visibility','hour_bef_pm10',target]
    else:
      use_columns = ['hour', 'hour_bef_temperature', 'hour_bef_precipitation','hour_bef_windspeed', 'hour_bef_humidity', 'hour_bef_visibility',target]

    #feature, target data split
    without_isna = pd.concat([train,test],axis=0).loc[:,use_columns].dropna()
    X_feature = without_isna.iloc[:,:-1]
    Y_train = without_isna[target]

    X_train, X_test, Y_train, Y_test = train_test_split(X_feature,Y_train, test_size=0.3, random_state=2021)

    # 원본 train data 안에 nan 값에 해당하는 Feature 추출 
    ozone_isna_train_df= train[train[target].isna()]

    # 원본 test data 안에 nan 값에 해당하는 Feature 추출
    ozone_isna_test_df= test[test[target].isna()]

    models = [
        ('ridge', lm.Ridge()),
        ('lasso', lm.Lasso()),
        ('elastic', lm.ElasticNet()),
        ('LassoLars', lm.LassoLars()),
        ('SGDRegressor', lm.SGDRegressor()),
        ('knn', KNeighborsRegressor(n_jobs = -1)),
    ]

    params = {
        'ridge': {
            'alpha': [0.01, 0.1, 1.0, 10, 100],
            'fit_intercept': [True, False],
            'normalize': [True, False],
        },
        'lasso': {
            'alpha': [0.1, 1.0, 10],
            'fit_intercept': [True, False],
            'normalize': [True, False],
        },
        'elastic': {
            'alpha': [0.1, 1.0, 10],
            'normalize': [True, False],
            'fit_intercept': [True, False],
        },
            'LassoLars': {
            'alpha': [0.1, 1.0, 10],
            'normalize': [True, False],
            'fit_intercept': [True, False],
        },
        'SGDRegressor': {
            'penalty': ['l1', 'l2'],
            'alpha': [0.001, 0.01, 0.1, 1.0, 10, 100],
            'fit_intercept': [True, False],
        },
        'knn': {
            "n_neighbors": range(2,7),
        }
    }
    
    #target 별 상관성이 상대적으로 낮은 변수들 제외하기 위한 x_colname 정의
    hour_bef_ozone_x_colname = [col for col in X_train_scaled.columns if col != 'hour_bef_precipitation' and col != 'hour_bef_ozone' and col != 'hour_bef_pm10' and col != 'hour_bef_pm2.5']
    hour_bef_pm10_x_colname =['hour_bef_humidity','hour_bef_visibility']
    hour_bef_pm25_x_colname = ['hour_bef_humidity','hour_bef_visibility','hour_bef_pm10'] 

    if target =="hour_bef_ozone":
      x_colname = hour_bef_ozone_x_colname
    elif target =="hour_bef_pm10":
      x_colname = hour_bef_pm10_x_colname
    else:
      x_colname = hour_bef_pm25_x_colname

    best_model, best_score = None, float('inf')
    for model_name, model in models:
        param_grid = params[model_name]
        grid = GridSearchCV(model, cv=5, n_jobs=-1, param_grid=param_grid)
        grid = grid.fit(X_train[x_colname], Y_train)

        model = grid.best_estimator_
        predictions = model.predict(X_test[x_colname])

        evaluation = evaluate(Y_test, predictions)
        score = evaluation['mse'][0]
        print(model_name, score)

        if score < best_score:
            best_score=score
            best_model = model

    print("최적의 모델:",best_model)
    print("최종 mse score:",best_score,"\n")
    best_model = grid.best_estimator_
    best_train_pred = best_model.predict(ozone_isna_train_df[x_colname])
    best_test_pred = best_model.predict(ozone_isna_test_df[x_colname])

    # 예측값을 train,test에 적용
    ozone_nan_train_index = train[train[target].isna()][target].index
    ozone_nan_test_index = test[test[target].isna()][target].index

    for pred,index in zip(best_train_pred,ozone_nan_train_index):
      train[target].fillna({index:pred}, inplace=True)
    for pred,index in zip(best_test_pred,ozone_nan_test_index):
      test[target].fillna({index:pred}, inplace=True)

targets =['hour_bef_ozone', 'hour_bef_pm10', 'hour_bef_pm2.5']
Nan_prediction(targets)

 

ridge 0.00022735981284121033
lasso 0.0003577431018933861
elastic 0.00032536467486597313
LassoLars 0.000357741115357743
SGDRegressor 0.0019727722132471728
knn 0.0002798427571351642
최적의 모델: Ridge(alpha=0.01, copy_X=True, fit_intercept=True, max_iter=None,
      normalize=True, random_state=None, solver='auto', tol=0.001)
최종 mse score: 0.00022735981284121033 

ridge 479.574083827637
lasso 480.9698920519781
elastic 480.99426051508937
LassoLars 480.97117933041903
SGDRegressor 3.519192783192444e+29
knn 514.6375284552845
최적의 모델: Ridge(alpha=0.01, copy_X=True, fit_intercept=True, max_iter=None,
      normalize=True, random_state=None, solver='auto', tol=0.001)
최종 mse score: 479.574083827637 

ridge 84.67710204218932
lasso 84.72469340584847
elastic 84.63747722153994
LassoLars 84.7245776238065
SGDRegressor 1.022157379699085e+30
knn 45.39878500823723
최적의 모델: KNeighborsRegressor(algorithm='auto', leaf_size=30, metric='minkowski',
                    metric_params=None, n_jobs=-1, n_neighbors=4, p=2,
                    weights='uniform')
최종 mse score: 45.39878500823723

각 target별 연관성이 있다고 판단한 변수들만 뽑아서 예측해 보았습니다. ozone변수를 보면 여러개의 X feature들과 연관성을 가지고 있습니다. 반면 pm10, pm 2.5는 영향력있는 변수가 적었기 때문에 예측성능이 좋지 못했습니다. 따라서 

ozone은 예측값을 사용하고 그외 pm10, pm 2.5은 예측값을 넣되 사용하지 하지않기로 결정했습니다. 왜냐하면 이 둘의 결측값의 수만 해도 100개 이상이며 최종 target인 count와 연관성이 적기 때문입니다. (아래 이미지를 참고해주세요)

 

 

 

* 다음 시간 예고

다음 시간부터는 count의 출퇴근 수요를 설명해줄 변수를 생성하려고 합니다. 

아래는 X feature들의 분포입니다.

colum=['hour_bef_temperature','hour_bef_windspeed','hour_bef_humidity','hour_bef_ozone','hour_bef_pm10','hour_bef_pm2.5']
for col in colum:
  train.groupby(by='hour').mean()[col].plot(label=col)
  plt.legend()
  plt.show()

 

 

*Count 분포

train.groupby(by='hour').mean()['count'].plot()

출퇴근 시간에 수요가 급등하는 모습을 보이고 있습니다. 이러한 특징적인 현상을 설명해줄 변수가 필요합니다. X feature들의 분포를 보면 몇몇 변수에서 퇴근시간 증가하는 모습을 보입니다. 그렇지만 증가하는 분포? 폭이 넓기 때문에 기대이상으로 퇴근 시간대의 수요를 설명해 주기란 어려울 것으로 생각했습니다. 특히 출근 시간대의 수요증가는 hour뿐 이기 때문에 새로운 변수의 필요성을 느끼게 됩니다. 아직은 이렇게 해서 만들어야 겠다! 라는 아이디어는 떠오르지 않아서 시간이 좀 걸리지 않을까 생각해 봅니다 

 

좀 더 개선된 내용을 가지고 다음 시간에 포스팅 하겠습니다. 감사합니다