Анализ и прогноза за цените на жилищния пазар с помощта на Cross Validation и Grid Search в няколко регресионни модела

В тази статия анализирам факторите, свързани с цените на жилищата в Мелбърн, и правя прогнози за цените на жилищата, използвайки няколко техники за машинно обучение: Линейна регресия, Регресия на Ридж, K-най-близки съседи (по-нататък KNN) и Дърво на решения. Използвайки методите на кръстосаното валидиране и техниките за търсене в решетка, намирам оптималните стойности за хиперпараметри във всеки модел и сравнявам резултатите, за да намеря най-добрия модел за машинно обучение за прогнозиране на цените на жилищата в Мелбърн.

Целият код за този проект и данните са тук.

Наборът от данни

Данните за този анализ са Пазарът на жилища в Мелбърн от набора от данни Kaggle. Общият брой на редовете и колоните е съответно 34 857 и 21. Колоните са както следва:

df = pd.read_csv('...\Melbourne_housing_FULL.csv')
df.columns.to_list()
['Suburb','Address','Rooms','Type','Method','SellerG','Date','Distance','Postcode','Bedroom2','Bathroom','Car','Landsize','BuildingArea','YearBuilt','CouncilArea','Latitude','Longitude','Regionname','Propertycount','Price']

Предварителна обработка на данни

В този раздел е обяснено накратко как се третират липсващите стойности и отклоненията.

Липсващи стойности

Използвайки библиотека missingno в Python, проверете липсващите стойности в набора от данни.

import missingno as msno
msno.bar(df)

  • Според коефициентите на корелация, Rooms е добър заместител на Bathroom и Car. След създаване на категоричен признак, указващ дали нивото на Rooms е високо, средно или нискоза всяка къща, медианите на Bathroom и Car са изчислено в рамките на група, която има същото ниво от Rooms. След това изчислените медиани се поставят за липсващите стойности в Bathroom и Car.
  • Редовете с липсващи стойности в Price се премахват.
  • Останалите функции с липсващи стойности се премахват.

Извънредни стойности

За извънредни стойности използвам метода Интерквартилен диапазон (IQR). Тази техника намира точки от данни, които попадат извън 1,5 пъти от интерквартилен диапазон над 3-тия квартил (Q3) и под 1-ви квартил (Q1), и премахва тези записи от анализа.

Q1 и Q3 за всяка цифрова характеристика и целта Price са както следва:

  • Price — $635 000 и $1 295 000
  • Rooms — 2 стаи и 4 стаи
  • Distance — 6,4 километра и 14 километра
  • Bathroom — 1 и 2 стаи.
  • Car — 1 и 2 петна

За Price точката за данни за 1,5 пъти IQR над Q3 (горен мустак) е $2 285 000, а точката за данни за 1,5 пъти IQR под Q1 (долна мустак) е -355 000 $. Разпределението на цената след премахване на отклоненията е както следва:

import seaborn as sns
import matplotlib.pyplot as plt
Q1 = df['Price'].quantile(0.25)
Q3 = df['Price'].quantile(0.75)
IQR = Q3-Q1
Lower_Whisker = Q1 - 1.5*IQR
Upper_Whisker = Q3 + 1.5*IQR
df = df[(df['Price']>Lower_Whisker)&(df['Price']<Upper_Whisker)]
plt.figure(figsize=(10,5))
sns.distplot(df['Price'],hist=True, kde=False, color='blue')
plt.ylabel('Counts')

Проучвателен анализ на данни (EDA)

Числени характеристики

В нашия анализ имаме четири числови характеристики: Rooms, Bathroom, Car и Distance. Най-добрият начин да видите връзките между числените характеристики и целта с един поглед е да начертаете диаграми на разсейване. Диаграма hexbinразделя областите на няколко hexbinв графиката и цветът на всеки hexbinозначава броя на точките с данни. Тъй като цветът на хексбинстава по-тъмен, има повече точки с данни в хексбин.

Нека да видим диаграмите на hexbinза разсейване, за да покажем връзките на числовите характеристики с Price.

import matplotlib.image as mpimg
JG1 = sns.jointplot('Rooms', 'Price', data=df, kind='hex', color='g')
JG2 = sns.jointplot('Bathroom', 'Price', data=df, kind='hex', color='b')
JG3 = sns.jointplot('Car', 'Price', data=df, kind='hex', color='r')
JG4 = sns.jointplot('Distance', 'Price', data=df, kind='hex', color='orange')
JG1.savefig('JG1.png')
plt.close(JG1.fig)
JG2.savefig('JG2.png')
plt.close(JG2.fig)
JG3.savefig('JG3.png')
plt.close(JG3.fig)
JG4.savefig('JG4.png')
plt.close(JG4.fig)
f, ax = plt.subplots(2,2,figsize=(20,16))
ax[0,0].imshow(mpimg.imread('JG1.png'))
ax[0,1].imshow(mpimg.imread('JG2.png'))
ax[1,0].imshow(mpimg.imread('JG3.png'))
ax[1,1].imshow(mpimg.imread('JG4.png'))
[ax.set_axis_off() for ax in ax.ravel()]
plt.tight_layout()

Графиките по-горе ясно показват положителните връзки на Rooms, Bathroom, Car с Price и отрицателната връзка на Distance с Price.

Категорични характеристики

Категориалните характеристики са Regionname и Type.

Име на регион

Regionname има 7 уникални стойности: Northern Metropolitan, Southern Metropolitan, Western Metropolitan, Eastern Metropolitan, South-Eastern Metropolitan, Northern Victoria, и Eизточна Виктория. Тези имена на региони са Избирателни региони на Виктория. По принцип тези региони са разделени на осем региона. В нашите данни няма данни за жилища, разположени в Западна Виктория.

Нека да видим графиката между Regionname и Price.

plt.figure(figsize=(12,6))
sns.boxplot('Regionname', 'Price', data=df, width=0.3, palette="Set2")
plt.xticks(rotation=45)
df['Regionname'].value_counts()

За да използвам тази функция в анализа, създавам манекени за тази функция и ги обединявам в набора от данни.

regionname = pd.get_dummies(df['Regionname'],drop_first=True)
df = pd.merge(df, regionname, left_index=True, right_index=True)
df.drop('Regionname', axis=1, inplace=True)

Тип

Нашият набор от данни разделя жилищните типове на три категории:

  • h — къща, вила, вила, полуфабрикат и тераса
  • u — единица и дуплекс
  • t — градска къща
plt.figure(figsize=(10,5))
sns.boxplot('Type', 'Price', data=df, width=0.3, palette="Set2")
df['Type'].value_counts()

Около 60% от наблюденията е тип h, а около 25% е тип u. В диаграмата типът h изглежда има най-голямата вариация. Най-скъпите и най-евтините къщи са от типа h.

Нека създадем манекени и за Type и да ги комбинираме в набора от данни.

house_type = pd.get_dummies(df['Type'], drop_first=True)
df = pd.merge(df,house_type, left_index=True, right_index=True)
df.drop('Type', axis=1, inplace=True)

Предсказуемо моделиране

Регресионните модели, използвани в този анализ, са линейна регресия, гребенна регресия, K- Най-близки съседи и Дърво на решения.

Основни прогнози

За да тествам ефективността на прогнозите за всеки регресионен модел, първо разделям данните на набор за обучение и тестване. Напасвам моделите, като първо използвам набора за обучение, и след това прогнозирам цените на жилищата, използвайки набора за тестване. Данните от тестването са определени да бъдат 30% от всички данни.

from sklearn.model_selection import train_test_split
from sklearn import metrics
from sklearn.model_selection import cross_validate
X=df.drop('Price', axis=1)
y=df['Price']
train_X, test_X, train_y, test_y = train_test_split(X,y,test_size=0.3, random_state=0)

За да измеря ефективността на прогнозите за всеки модел, използвам два показателя за ефективност: R²(Коефициент на определяне) и MSE (Средна квадратна грешка)

Първият е коефициентът на определяне, който обикновено се изразява като R². Коефициентът на определяне е съотношението на дисперсията на дадена цел, обяснена или прогнозирана от модел, към общата дисперсия на целта. Тя варира от 0 до 1 и тъй като стойността е по-близо до 1, моделът обяснява или прогнозира дисперсията на целта по-добре.

  • MSE

Вторият показател е средна квадратична грешка (MSE). MSE е средната стойност на квадратната разлика между оценените или прогнозираните стойности и действителните стойности на дадена цел. Това винаги е по-голямо от нула. По-ниската стойност на MSE показва по-висока точност на прогнозите на модела. В този анализ използвам корен квадратен от този показател (RMSE).

За удобство нека създадем функция за изчисляване на R² и RMSE и функция за сравняване на разпределенията на действителните стойности и прогнозираните стойности за всеки модел.

def Predictive_Model(estimator):
    estimator.fit(train_X, train_y)
    prediction = estimator.predict(test_X)
    print('R_squared:', metrics.r2_score(test_y, prediction))
    print('Square Root of MSE:',np.sqrt(metrics.mean_squared_error(test_y, prediction)))
    plt.figure(figsize=(10,5))
    sns.distplot(test_y, hist=True, kde=False)
    sns.distplot(prediction, hist=True, kde=False)
    plt.legend(labels=['Actual Values of Price', 'Predicted Values of Price'])
    plt.xlim(0,)
def FeatureBar(model_Features, Title, yLabel):
    plt.figure(figsize=(10,5))
    plt.bar(df.columns[df.columns!='Price'].values, model_Features)
    plt.xticks(rotation=45)
    plt.title(Title)
    plt.ylabel(yLabel)

Линейна регресия

from sklearn.linear_model import LinearRegression
lr = LinearRegression()
Predictive_Model(lr)

R² на линейната регресия е 0,6146, което означава, че около 60% от дисперсията на цените на жилищата в данните може да бъде предвидена от модела. RMSE е 264 465. Това означава, че за всички прогнози за набора от тестове, средната разлика за всяка прогноза е $264 465.

Регресия на билото

from sklearn.linear_model import Ridge
rr = Ridge(alpha=100)
Predictive_Model(rr)

Резултатът по-горе се получава с параметър за регулиране (alpha), равен на 100. R² на този модел е 0,6133, а RMSE е 264 920.

K-най-близки съседи (KNN)

from sklearn.neighbors import KNeighborsRegressor
knn = KNeighborsRegressor(n_neighbors=5)
Predictive_Model(knn)

Резултатът по-горе се получава с броя на съседите, равен на 5. R² на този модел е 0,7053, а RMSE е 231,250.

Дърво на решения

from sklearn.tree import DecisionTreeRegressor
dt = DecisionTreeRegressor(max_depth=15, random_state=0)
Predictive_Model(dt)

В модела на дървото на решенията max depth е един от факторите за предотвратяване на проблема с прекомерното приспособяване на модела. Тъй като дълбочината на дървото е по-голяма, дървото има повече клони и става по-голямо. Тъй като дървото има повече клони, прогнозата за набора за обучение може да бъде по-точна. Въпреки това, има по-голяма вариация в прогнозирането на набора за тестване. Следователно оптималната настройка на max depth е важна, за да се избегне проблемът с пренапасването. В горния пример max depth е настроен на 15. R² на този модел е 0,6920, а RMSE е 236 424.

Обобщение на ефективността

regressor = ['Linear Regression', 'Ridge Regression', 'KNN', 'Decision Tree']
models = [LinearRegression(), Ridge(alpha=100), KNeighborsRegressor(n_neighbors=5), DecisionTreeRegressor(max_depth=15, random_state=0)]
R_squared = []
RMSE = []
for m in models:
    m.fit(train_X, train_y)
    prediction_m = m.predict(test_X)
    r2 = metrics.r2_score(test_y, prediction_m)
    rmse = np.sqrt(metrics.mean_squared_error(test_y, prediction_m))
    R_squared.append(r2)
    RMSE.append(rmse)
basic_result = pd.DataFrame({'R squared':R_squared,'RMSE':RMSE}, index=regressor)
basic_result

В таблицата по-горе KNN изглежда оптималният модел за прогнозиране на цените на жилищата в Мелбърн. Въпреки това е твърде рано да се правят заключения, тъй като има още неща, които трябва да бъдат разгледани. Първият е, че използваме само една конкретна селекция от комплекти за обучение и тестване, а вторият е, че за всеки модел избрахме една конкретна стойност за всеки хипер параметър. За да получим стабилен резултат, покриващ тези проблеми, трябва също да преминем през процеса на кръстосано валидиране и търсене в мрежата.

Кръстосано валидиране и търсене в мрежа

Кръстосано валидиране (CV) е процедура за повторно вземане на проби, когато количеството данни е ограничено. Това на случаен принцип разделя всички данни на K-сгъвания, напасва модел с помощта на (K-1) сгъвания, валидира модела с помощта на оставащото сгъване и след това оценява производителността чрез показатели. След това CV повтаря целия този процес, докато всяко K-кратно се използва като набор за тестване. Средната стойност на K-броя на резултатите на метрика е крайната оценка на ефективността за модела.

Търсене в мрежа е процес на настройка на хиперпараметри, за да се намерят оптималните стойности на параметрите за модел. Резултатите от прогнозата могат да варират в зависимост от конкретните стойности на параметрите. Техниката за търсене в мрежата прилага всички възможни кандидати за параметрите, за да намери оптималния, за да даде най-добрите прогнози за модела.

Линейна регресия

scoring={'R_squared':'r2','MSE':'neg_mean_squared_error'}
def CrossVal(estimator):
    scores = cross_validate(estimator, X, y, cv=10, scoring=scoring)
    r2 = scores['test_R_squared'].mean()
    mse = abs(scores['test_Square Root of MSE'].mean())
    print('R_squared:', r2)
    print('Square Root of MSE:', np.sqrt(mse))
CrossVal(LinearRegression())
R_squared: 0.5918115585795747
Square Root of MSE: 269131.0885647736

Тъй като линейната регресия няма хипер параметър в нашия анализ, тук се извършва само CV. Броят на гънките в автобиографията е зададен на 10. Средната стойност на R² е 0,5918, а RMSE е 269131.

Регресия на билото

Параметърът за регулиране в ръбовата регресия се изразява като alpha в sklearn. Тъй като GridSearchCV в sklearn включва процеса на кръстосано валидиране, процесът на изпълнение на cross_validate е пропуснат. Наборът от мрежата за alpha е настроен да бъде [0.01, 0.1, 1, 10, 100, 1000, 10000] тук.

from sklearn.model_selection import GridSearchCV
def GridSearch(estimator, Features, Target, param_grid):
    for key, value in scoring.items():
        grid = GridSearchCV(estimator, param_grid, cv=10, scoring=value)
        grid.fit(Features,Target)
        print(key)
        print('The Best Parameter:', grid.best_params_)
        if grid.best_score_ > 0:
            print('The Score:', grid.best_score_)
        else:
            print('The Score:', np.sqrt(abs(grid.best_score_)))
        print()
param_grid = {'alpha':[0.01, 0.1, 1, 10, 100, 1000, 10000]}
GridSearch(Ridge(), X, y, param_grid)
R_squared
The Best Parameter: {'alpha': 10}
The Score: 0.5918404945235951
Square Root of MSE
The Best Parameter: {'alpha': 10}
The Score: 269125.20208461734

Резултатът показва, че най-добрата стойност на alpha е 10. R² и RMSE са съответно 0,5918 и 269125 под alpha = 10.

K-най-близки съседи

Хиперпараметърът за KNN, който използваме в този анализ, е броят на най-близките съседи (n_neighbors). Диапазонът за мрежата е целите числа от 5 до 25.

param_grid = dict(n_neighbors=np.arange(5,26))
GridSearch(KNeighborsRegressor(), X, y, param_grid)
R_squared
The Best Parameter: {'n_neighbors': 16}
The Score: 0.6973921821195777
Square Root of MSE
The Best Parameter: {'n_neighbors': 16}
The Score: 232900.0204190322

Оптималният брой на n_neighbors е 16. R² е 0,6974, а RMSE е 232900. Можем да видим как 16 е оптималната стойност за n_neighbors в нашия анализ, като разгледаме кривата на валидиране.

from sklearn.model_selection import validation_curve
def ValidationCurve(estimator, Features, Target, param_name, Name_of_HyperParameter, param_range):
    
    train_score, test_score = validation_curve(estimator, Features, Target, param_name, param_range,cv=10,scoring='r2')
    Rsqaured_train = train_score.mean(axis=1)
    Rsquared_test= test_score.mean(axis=1)
    
    plt.figure(figsize=(10,5))
    plt.plot(param_range, Rsqaured_train, color='r', linestyle='-', marker='o', label='Training Set')
    plt.plot(param_range, Rsquared_test, color='b', linestyle='-', marker='x', label='Testing Set')
    plt.legend(labels=['Training Set', 'Testing Set'])
    plt.xlabel(Name_of_HyperParameter)
    plt.ylabel('R_squared')
ValidationCurve(KNeighborsRegressor(), X, y, 'n_neighbors', 'K-Neighbors',np.arange(5,26))

Дърво на решения

В модела на дървото на решенията може да има няколко хиперпараметъра, които трябва да се вземат предвид. В нашия анализ само max_depth е опцията за хиперпараметъра. Обхватът на max_depth за проверка е целите числа от 2 до 14.

param_grid=dict(max_depth=np.arange(2,15))
GridSearch(DecisionTreeRegressor(), X, y, param_grid)
R_squared
The Best Parameter: {'max_depth': 9}
The Score: 0.6844562874572124
Square Root of MSE
The Best Parameter: {'max_depth': 9}
The Score: 237708.76352194021

Резултатът показва, че оптималната стойност за max_depth е 9. R² е 0,6845, а RMSE е 237708 под max_depth=9. Това се потвърждава и от кривата на валидиране.

ValidationCurve(DecisionTreeRegressor(), X, y, 'max_depth', 'Maximum Depth', np.arange(4,15))

Обобщение на кръстосаното валидиране

Таблицата и графиките по-долу показват резултатите на R² за всеки кръг на тестване в CV. Тъй като cv е настроено да бъде 10, имаме 10 кръга тестове.

lr_scores = cross_validate(LinearRegression(), X, y, cv=10, scoring='r2')
rr_scores = cross_validate(Ridge(alpha=10), X, y, cv=10, scoring='r2')
knn_scores = cross_validate(KNeighborsRegressor(n_neighbors=16), X, y, cv=10, scoring='r2')
dt_scores = cross_validate(DecisionTreeRegressor(max_depth=9, random_state=0), X, y, cv=10, scoring='r2')
lr_test_score = lr_scores.get('test_score')
rr_test_score = rr_scores.get('test_score')
knn_test_score = knn_scores.get('test_score')
dt_test_score = dt_scores.get('test_score')
box= pd.DataFrame({'Linear Regression':lr_test_score, 'Ridge Regression':rr_test_score, 'K-Nearest Neighbors':knn_test_score, 'Decision Tree':dt_test_score})
box.index = box.index + 1
box.loc['Mean'] = box.mean()
box

Според резултата в таблицата най-добрият модел за машинно обучение в нашия анализ е KNN, тъй като средната стойност на резултатите за всеки кръг е най-високата за KNN.

f,ax=plt.subplots(1,2, figsize=(12,5))
sns.boxplot(data=box.drop(box.tail(1).index), width=0.3, palette="Set2", ax=ax[0])
ax[0].set_ylabel('R squared')
sns.lineplot(data=box.drop(box.tail(1).index), palette="Set2", ax=ax[1])
ax[1].set_xticks(np.arange(1,11,1))
ax[1].set_xlabel('K-th Fold')

Каретата и линейните диаграми по-горе показват разпределенията и промените на резултатите за всеки модел. Моделът на дървото на решенията, както и KNN показват добро представяне в нашия анализ. Линейната и ръбовата регресия не показват значителна разлика в тяхното представяне.