哈喽，大家好。
今天給大家介紹一個非常適合新手入門的機器學習實戰案例。
這是一個房價預測的案例，來源於 Kaggle 網站，是很多算法初學者的第一道競賽題目。
該案例有著解機器學習問題的完整流程，包含EDA、特征工程、模型訓練、模型融合等。
房價預測流程
下面跟著我，來學習一下該案例。
沒有啰嗦的文字，沒有多余的代碼，只有通俗的講解。
探索性數據分析(Exploratory Data Analysis，簡稱EDA) 的目的是讓我們對數據集有充分的了解。在這一步，我們探索的內容如下：
EDA內容
1.1 輸入數據集
train = pd.read_csv('./data/train.csv')
test = pd.read_csv('./data/test.csv')
訓練樣本
train和test分別是訓練集和測試集，分別有 1460 個樣本，80 個特征。
SalePrice列代表房價，是我們要預測的。
1.2 房價分布
因為我們任務是預測房價，所以在數據集中核心要關注的就是房價（SalePrice) 一列的取值分布。
sns.distplot(train['SalePrice']);
房價取值分布
從圖上可以看出，SalePrice列峰值比較陡，並且峰值向左偏。
也可以直接調用skew()和kurt()函數計算SalePrice具體的偏度和峰度值。
對於偏度和峰度都比較大的情況，建議對SalePrice列取log()進行平滑。
1.3 與房價相關的特征
了解完SalePrice的分布後，我們可以計算 80 個特征與SalePrice的相關關系。
重點關注與SalePrice相關性最強的 10 個特征。
# 計算列之間相關性
corrmat = train.corr()
# 取 top10
k = 10
cols = corrmat.nlargest(k, 'SalePrice')['SalePrice'].index
# 繪圖
cm = np.corrcoef(train[cols].values.T)
sns.set(font_scale=1.25)
hm = sns.heatmap(cm, cbar=True, annot=True, square=True, fmt='.2f', annot_kws={'size': 10}, yticklabels=cols.values, xticklabels=cols.values)
plt.show()
與SalePrice高度相關的特征
OverallQual（房子材料和裝飾）、GrLivArea（地上居住面積）、GarageCars（車庫容量）和 TotalBsmtSF（地下室面積）跟SalePrice有很強的相關性。
這些特征在後面做特征工程時也會重點關注。
1.4 剔除離群樣本
由於數據集樣本量很少，離群點不利於我們後面訓練模型。
所以需要計算每個數值特性的離群點，剔除掉離群次數最多的樣本。
# 獲取數值型特征
numeric_features = train.dtypes[train.dtypes != 'object'].index
# 計算每個特征的離群樣本
for feature in numeric_features:
outs = detect_outliers(train[feature], train['SalePrice'],top=5, plot=False)
all_outliers.extend(outs)
# 輸出離群次數最多的樣本
print(Counter(all_outliers).most_common())
# 剔除離群樣本
train = train.drop(train.index[outliers])
detect_outliers()是自定義函數，用sklearn庫的LocalOutlierFactor算法計算離群點。
到這裡， EDA 就完成了。最後，將訓練集和測試集合並，進行下面的特征工程。
y = train.SalePrice.reset_index(drop=True)
train_features = train.drop(['SalePrice'], axis=1)
test_features = test
features = pd.concat([train_features, test_features]).reset_index(drop=True)
features合並了訓練集和測試集的特征，是我們下面要處理的數據。
特征工程
2.1 校正特征類型
MSSubClass（房屋類型）、YrSold（銷售年份）和MoSold（銷售月份）是類別型特征，只不過用數字來表示，需要將它們轉成文本特征。
features['MSSubClass'] = features['MSSubClass'].apply(str)
features['YrSold'] = features['YrSold'].astype(str)
features['MoSold'] = features['MoSold'].astype(str)
2.2 填充特征缺失值
填充缺失值沒有統一的標准，需要根據不同的特征來決定按照什麼樣的方式來填充。
# Functional：文檔提供了典型值 Typ
features['Functional'] = features['Functional'].fillna('Typ') #Typ 是典型值
# 分組填充需要按照相似的特征分組，取眾數或中位數
# MSZoning（房屋區域）按照 MSSubClass（房屋）類型分組填充眾數
features['MSZoning'] = features.groupby('MSSubClass')['MSZoning'].transform(lambda x: x.fillna(x.mode()[0]))
#LotFrontage（到接到舉例）按Neighborhood分組填充中位數
features['LotFrontage'] = features.groupby('Neighborhood')['LotFrontage'].transform(lambda x: x.fillna(x.median()))
# 車庫相關的數值型特征，空代表無，使用0填充空值。
for col in ('GarageYrBlt', 'GarageArea', 'GarageCars'):
features[col] = features[col].fillna(0)
2.3 偏度校正
跟探索SalePrice列類似，對偏度高的特征進行平滑。
# skew()方法，計算特征的偏度（skewness）。
skew_features = features[numeric_features].apply(lambda x: skew(x)).sort_values(ascending=False)
# 取偏度大於 0.15 的特征
high_skew = skew_features[skew_features > 0.15]
skew_index = high_skew.index
# 處理高偏度特征，將其轉化為正態分布，也可以使用簡單的log變換
for i in skew_index:
features[i] = boxcox1p(features[i], boxcox_normmax(features[i] + 1))
2.4 特征刪除和新增
對於幾乎都是缺失值，或單一取值占比高（99.94%）的特征可以直接刪除。
features = features.drop(['Utilities', 'Street', 'PoolQC',], axis=1)
同時，可以融合多個特征，生成新特征。
有時候模型很難學習到特征之間的關系，手動融合特征可以降低模型學習難度，提升效果。
# 將原施工日期和改造日期融合
features['YrBltAndRemod']=features['YearBuilt']+features['YearRemodAdd']
# 將地下室面積、1樓、2樓面積融合
features['TotalSF']=features['TotalBsmtSF'] + features['1stFlrSF'] + features['2ndFlrSF']
可以發現，我們融合的特征都是與SalePrice強相關的特征。
最後簡化特征，對分布單調的特征（如：100個數據中有99個的數值是0.9，另1個是0.1），進行01處理。
features['haspool'] = features['PoolArea'].apply(lambda x: 1 if x > 0 else 0)
features['has2ndfloor'] = features['2ndFlrSF'].apply(lambda x: 1 if x > 0 else 0)
2.6 生成最終訓練數據
到這裡特征工程就做完了， 我們需要從features中將訓練集和測試集重新分離出來，構造最終的訓練數據。
X = features.iloc[:len(y), :]
X_sub = features.iloc[len(y):, :]
X = np.array(X.copy())
y = np.array(y)
X_sub = np.array(X_sub.copy())
因為SalePrice是數值型且是連續的，所以需要訓練一個回歸模型。
3.1 單一模型
首先以嶺回歸（Ridge） 為例，構造一個k折交叉驗證模型。
from sklearn.linear_model import RidgeCV
from sklearn.pipeline import make_pipeline
from sklearn.model_selection import KFold
kfolds = KFold(n_splits=10, shuffle=True, random_state=42)
alphas_alt = [14.5, 14.6, 14.7, 14.8, 14.9, 15, 15.1, 15.2, 15.3, 15.4, 15.5]
ridge = make_pipeline(RobustScaler(), RidgeCV(alphas=alphas_alt, cv=kfolds))
嶺回歸模型有一個超參數alpha，而RidgeCV的參數名是alphas，代表輸入一個超參數alpha數組。在擬合模型時，會從alpha數組中選擇表現較好某個取值。
由於現在只有一個模型，無法確定嶺回歸是不是最佳模型。所以我們可以找一些出場率高的模型多試試。
# lasso
lasso = make_pipeline(
RobustScaler(),
LassoCV(max_iter=1e7, alphas=alphas2, random_state=42, cv=kfolds))
#elastic net
elasticnet = make_pipeline(
RobustScaler(),
ElasticNetCV(max_iter=1e7, alphas=e_alphas, cv=kfolds, l1_ratio=e_l1ratio))
#svm
svr = make_pipeline(RobustScaler(), SVR(
C=20,
epsilon=0.008,
gamma=0.0003,
))
#GradientBoosting（展開到一階導數）
gbr = GradientBoostingRegressor(...)
#lightgbm
lightgbm = LGBMRegressor(...)
#xgboost（展開到二階導數）
xgboost = XGBRegressor(...)
有了多個模型，我們可以再定義一個得分函數，對模型評分。
#模型評分函數
def cv_rmse(model, X=X):
rmse = np.sqrt(-cross_val_score(model, X, y, scoring="neg_mean_squared_error", cv=kfolds))
return (rmse)
以嶺回歸為例，計算模型得分。
score = cv_rmse(ridge)
print("Ridge score: {:.4f} ({:.4f})\n".format(score.mean(), score.std()), datetime.now(), ) #0.1024
運行其他模型發現得分都差不多。
這時候我們可以任選一個模型，擬合，預測，提交訓練結果。還是以嶺回歸為例
# 訓練模型
ridge.fit(X, y)
# 模型預測
submission.iloc[:,1] = np.floor(np.expm1(ridge.predict(X_sub)))
# 輸出測試結果
submission = pd.read_csv("./data/sample_submission.csv")
submission.to_csv("submission_single.csv", index=False)
submission_single.csv是嶺回歸預測的房價，我們可以把這個結果上傳到 Kaggle 網站查看結果的得分和排名。
3.2 模型融合-stacking
有時候為了發揮多個模型的作用，我們會將多個模型融合，這種方式又被稱為集成學習。
stacking 是一種常見的集成學習方法。簡單來說，它會定義個元模型，其他模型的輸出作為元模型的輸入特征，元模型的輸出將作為最終的預測結果。
stacking
這裡，我們用mlextend庫中的StackingCVRegressor模塊，對模型做stacking。
stack_gen =
StackingCVRegressor(
regressors=(ridge, lasso, elasticnet, gbr, xgboost, lightgbm),
meta_regressor=xgboost,
use_features_in_secondary=True)
訓練、預測的過程與上面一樣，這裡不再贅述。
3.3 模型融合-線性融合
多模型線性融合的思想很簡單，給每個模型分配一個權重（權重加和=1），最終的預測結果取各模型的加權平均值。
# 訓練單個模型
ridge_model_full_data = ridge.fit(X, y)
lasso_model_full_data = lasso.fit(X, y)
elastic_model_full_data = elasticnet.fit(X, y)
gbr_model_full_data = gbr.fit(X, y)
xgb_model_full_data = xgboost.fit(X, y)
lgb_model_full_data = lightgbm.fit(X, y)
svr_model_full_data = svr.fit(X, y)
models = [
ridge_model_full_data, lasso_model_full_data, elastic_model_full_data,
gbr_model_full_data, xgb_model_full_data, lgb_model_full_data,
svr_model_full_data, stack_gen_model
]
# 分配模型權重
public_coefs = [0.1, 0.1, 0.1, 0.1, 0.15, 0.1, 0.1, 0.25]
# 線性融合，取加權平均
def linear_blend_models_predict(data_x,models,coefs, bias):
tmp=[model.predict(data_x) for model in models]
tmp = [c*d for c,d in zip(coefs,tmp)]
pres=np.array(tmp).swapaxes(0,1)
pres=np.sum(pres,axis=1)
return pres
到這裡，房價預測的案例我們就講解完了，大家可以自己運行一下，看看不同方式訓練出來的模型效果。
回顧整個案例會發現，我們在數據預處理和特征工程上花費了很大心思，雖然機器學習問題模型原理比較難學，但實際過程中往往特征工程花費的心思最多。
完整源代碼關注公眾號：Python源碼 即可獲取
