比赛是Kaggle举办的，要求选手根据爱荷华州埃姆斯市住宅的各种信息来预测房屋的最终价格。我认为这是一个回归问题，从数据清洗中的对异常值和缺失值的处理到特征工程中对特征的选择和设计，再到建模预测中对回归算法的选择都影响到结果的准确性。
比赛数据包括2个数据集，分别是训练集和测试集。其中训练集中有80个字段，测试集中有79个字段，ID是每一个房屋信息的编号，训练集比测试集多一列SalePrice的数据，而这个销售价格就是我们预测的目标。各主要字段的含义如下所示： SalePrice-物业的销售价格（美元） MSSubClass：建筑类 MSZoning：常规分区分类 LotFrontage：连接到物业的街道的线性英尺 LotArea：平方英尺大小 街道：道路通行类型 胡同：胡同通道的类型 LotShape：属性的一般形状 LandContour：物业的平坦度 实用程序：可用的实用程序类型 LotConfig：批次配置 LandSlope：物业的坡度 邻里：艾姆斯市区范围内的地理位置 条件1：接近主干道或铁路 条件2：接近主要道路或铁路（如果有） BldgType：住宅类型 HouseStyle：住宅风格 OverallQual：总体材料和加工质量 OverallCond：总体状况的评价 YearBuilt：原始施工日期 YearRemodAdd：改型日期 RoofStyle：屋顶类型 RoofMatl：屋顶材料 Exterior1st：房屋外墙 Exterior2nd：房屋的外墙覆盖物（如果一种以上的材料） MasVnrType：石工饰面类型 MasVnrArea：砌面贴面面积（平方英尺） ExterQual：外部材料质量 ExterCond：外部材料的当前状态 基金会：基金会的类型 BsmtQual：地下室的高度 BsmtCond：地下室的一般状况 BsmtExposure：罢工或花园水平的地下室墙 BsmtFinType1：地下室成品区域的质量 BsmtFinSF1：1型成品平方英尺 BsmtFinType2：第二个完成区域的质量（如果存在） BsmtFinSF2：2型成品平方英尺 BsmtUnfSF：地下室面积未完成的平方英尺 TotalBsmtSF：地下室总面积 加热方式：加热方式 加热质量控制：加热质量和条件 CentralAir：中央空调 电气：电气系统 1stFlrSF：第一层平方英尺 2ndFlrSF：二楼平方英尺 LowQualFinSF：低质量成品平方英尺（所有楼层） GrLivArea：地面（地面）以上的居住面积平方英尺 BsmtFullBath：地下室全浴室 BsmtHalfBath：地下室半浴室 FullBath：满级浴室 HalfBath：地上半浴 卧室：地下室以上的卧室数量 厨房：厨房数量 KitchenQual：厨房质量 TotRmsAbvGrd：上等客房总数（不包括浴室） 功能性：家庭功能性等级 壁炉：壁炉数量 FireplaceQu：壁炉质量 车库类型：车库位置 GarageYrBlt：年车库已建成 GarageFinish：车库的内部装饰 GarageCars：车库中车库的大小 GarageArea：车库的大小（平方英尺） GarageQual：车库质量 GarageCond：车库条件 PavedDrive：铺砌的车道 WoodDeckSF：木制甲板面积（平方英尺） OpenPorchSF：开放式阳台面积（平方英尺） EnclosedPorch：封闭的门廊面积（以平方英尺为单位） 3SsnPorch：三季门廊面积（以平方英尺为单位） ScreenPorch：屏幕门廊面积（以平方英尺为单位） PoolArea：平方英尺的游泳池面积 PoolQC：泳池质量 围栏：围栏质量 MiscFeature：杂项功能未包括在其他类别 MiscVal：杂项功能的$ Value MoSold：已售一个月 年销售额：已售年份 SaleType：销售类型 SaleCondition：销售条件
首先我们导入一些需要的工具包
import numpy as np # linear algebra import pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv) %matplotlib inline import matplotlib.pyplot as plt # Matlab-style plotting import seaborn as sns color = sns.color_palette() sns.set_style('darkgrid') import warnings def ignore_warn(*args, **kwargs): pass warnings.warn = ignore_warn #ignore annoying warning (from sklearn and seaborn) from scipy import stats from scipy.stats import norm, skew #for some statistics然后我们读取数据集
train = pd.read_csv('C:/Users/CSC/Desktop/train.csv') test = pd.read_csv('C:/Users/CSC/Desktop/test.csv')我们可以分别看一下训练集和测试集
train.head(5) test.head(5)首先我们检查不同特征的缺失情况。id这一列对于模型预测是没有用的，我们先将它在训练集和测试集中去掉，然后合并train test数据一起处理。
train_ID = train['Id'] test_ID = test['Id'] train.drop("Id", axis = 1, inplace = True) test.drop("Id", axis = 1, inplace = True) ntrain = train.shape[0] ntest = test.shape[0] y_train = train.SalePrice.values all_data = pd.concat((train, test)).reset_index(drop=True) all_data.drop(['SalePrice'], axis=1, inplace=True)接着我们计算每个字段的缺失率
all_data_na = (all_data.isnull().sum() / len(all_data)) * 100 all_data_na = all_data_na.drop(all_data_na[all_data_na == 0].index).sort_values(ascending=False)[:30] missing_data = pd.DataFrame({'Missing Ratio' :all_data_na}) missing_data.head(20)通过图表我们可以更直观地看到缺失率的情况
f, ax = plt.subplots(figsize=(15, 12)) sns.barplot(x=all_data_na.index, y=all_data_na) plt.xlabel('Features', fontsize=15) plt.ylabel('Percent of missing values', fontsize=15) plt.title('Percent missing data by feature', fontsize=15) plt.xticks(rotation='90')基于上面的检查，我们就探索完成了数据集中不同特征的缺失情况，对缺失值的处理，我们在后面的数据处理中进行，下面我们来探索不同特征和房价的相关性情况，同时对可能存在的异常值进行识别。
首先我们通过热力图来显示变量之间的相关性，然后选择和房价相关性最强的10个变量查看相关系数。
corrmat = train.corr() f, ax = plt.subplots(figsize=(12, 9)) sns.heatmap(corrmat, vmax=.8, square=True) 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()通过上面的热力图，我们可以有如下的结论：
1.‘OverallQual’, ‘GrLivArea’ 和 ‘TotalBsmtSF’与房价的相关性很强，可以后面再深入探索 2.一些变量之间也存在很高的相关性，比如’GarageCars’（车库能放多少量车） 和 ‘GarageArea’ （车库面积）这两个变量本身的相关性很强，这会涉及到多重共线性的问题，对于这两个变量，我们可以去掉一个，只留一个就可以了。还有’TotalBsmtSF’ 和 ‘1stFloor’ 似乎都是代表地下室面积，它们也有很强的相关性。后面我们在数据处理中会考虑进行处理。
接着我们查看一下相关性较强的几个变量的异常值，对于连续型变量我们采取散点图的方式，而离散型变量我们通过箱型图的方式来判断异常值。 （1）GrLivArea代表含义是居住面积，发现和房价有比较明显的正相关关系。我们可以在右下角看到两个明显的outlier，价格很低，但面积巨大，我们可以考虑删除掉这两个异常值。
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()（2）TotalBsmtSF含义为地下室面积，发现地下室面积与房价似乎有更强的潜在线性关系，同时在右侧似乎也有一个异常值存在。
var = 'TotalBsmtSF' data = pd.concat([train['SalePrice'], train[var]], axis=1) data.plot.scatter(x=var, y='SalePrice', ylim=(0,800000));（3）GarageArea含义为车库面积，可以看到车库面积和房价也存在一定的正相关关系。
var = 'GarageArea' data = pd.concat([train['SalePrice'], train[var]], axis=1) data.plot.scatter(x=var, y='SalePrice', ylim=(0,800000));接下来是离散型变量 （4）首先看一下房屋质量和房价之间的关系，可以看到随着房屋整体质量的变好，房屋的整体价格也在逐渐提高，是一个比较强力的预测变量。
var = 'OverallQual' data = pd.concat([train['SalePrice'], train[var]], axis=1) f, ax = plt.subplots(figsize=(8, 6)) fig = sns.boxplot(x=var, y="SalePrice", data=data) fig.axis(ymin=0, ymax=800000);（5）再看一下房屋的建造时间和价格的关系，可以看到房屋的建造时间虽然没有和房价有明显的线性关系，但最近建造的新房屋，整体的价格会比较高，这个在后面的特征工程中我们会采取虚拟变量即哑变量的形式来描述。
var = 'YearBuilt' data = pd.concat([train['SalePrice'], train[var]], axis=1) f, ax = plt.subplots(figsize=(16, 8)) fig = sns.boxplot(x=var, y="SalePrice", data=data) fig.axis(ymin=0, ymax=800000); plt.xticks(rotation=90);通过可视化探索，我们对数据可以建立起基本的了解，下一步我们对数据的质量做一些处理。首先是异常值，我们已经发现在房屋面积和地下室面积这两个特征里面可能存在异常值，可以把它们删掉。
对于房屋面积
train = train.drop(train[(train['GrLivArea']>4000) & (train['SalePrice']<300000)].index)对于地下室面积
train = train.drop(train[(train['TotalBsmtSF']>5000) & (train['SalePrice']<200000)].index)再来看一下异常值的分布情况
all_data_na = (all_data.isnull().sum() / len(all_data)) * 100 all_data_na = all_data_na.drop(all_data_na[all_data_na == 0].index).sort_values(ascending=False)[:30] missing_data = pd.DataFrame({'Missing Ratio' :all_data_na}) missing_data.head(50)PoolQC、MiscFeature、Alley的缺失值都在90%以上，可以考虑直接删掉这些特征
all_data= all_data.drop('PoolQC',axis=1) all_data = all_data.drop('MiscFeature',axis=1) all_data= all_data.drop('Alley',axis=1)Fence栅栏，FireplaceQu壁炉，如果数据缺失的话可能是代表房屋没有栅栏或壁炉， 我们用None填补缺失值，代表没有
all_data["Fence"] = all_data["Fence"].fillna("None") all_data["FireplaceQu"] = all_data["FireplaceQu"].fillna("None")LotFrontage代表房屋前街道的长度，从现实中考虑，房屋前街道的长度可能是会和房屋所在一个街区的房屋相同，所以我们可以考虑同一个街区房屋LotFrontage的均值来填补缺失值
all_data["LotFrontage"] = all_data.groupby("Neighborhood")["LotFrontage"].transform( lambda x: x.fillna(x.median()))Garage相关的车库变量，注意到这些变量的缺失比例是完全相同的，缺失这些变量的房屋可能是没有车库，用None填充
for col in ('GarageType', 'GarageFinish', 'GarageQual', 'GarageCond'): all_data[col] = all_data[col].fillna('None')同样猜测缺失值缺失的原因可能是因为房屋没有车库，对于连续型变量用0填充
for col in ('GarageYrBlt', 'GarageArea', 'GarageCars'): all_data[col] = all_data[col].fillna(0)地下室相关连续变量，缺失同样认为房屋可能是没有地下室，用0填充
for col in ('BsmtFinSF1', 'BsmtFinSF2', 'BsmtUnfSF','TotalBsmtSF', 'BsmtFullBath', 'BsmtHalfBath'): all_data[col] = all_data[col].fillna(0)地下室相关离散变量，同理用None填充
for col in ('BsmtQual', 'BsmtCond', 'BsmtExposure', 'BsmtFinType1', 'BsmtFinType2'): all_data[col] = all_data[col].fillna('None')Mas为砖石结构相关变量，缺失值我们同样认为是没有砖石结构，用0和none填补缺失值
all_data["MasVnrType"] = all_data["MasVnrType"].fillna("None") all_data["MasVnrArea"] = all_data["MasVnrArea"].fillna(0)对于其他缺失值处理的思路，和上面类似
all_data['MSZoning'] = all_data['MSZoning'].fillna(all_data['MSZoning'].mode()[0]) all_data = all_data.drop(['Utilities'], axis=1) all_data["Functional"] = all_data["Functional"].fillna("Typ") all_data['Electrical'] = all_data['Electrical'].fillna(all_data['Electrical'].mode()[0]) all_data['KitchenQual'] = all_data['KitchenQual'].fillna(all_data['KitchenQual'].mode()[0]) all_data['Exterior1st'] = all_data['Exterior1st'].fillna(all_data['Exterior1st'].mode()[0]) all_data['Exterior2nd'] = all_data['Exterior2nd'].fillna(all_data['Exterior2nd'].mode()[0]) all_data['SaleType'] = all_data['SaleType'].fillna(all_data['SaleType'].mode()[0]) all_data['MSSubClass'] = all_data['MSSubClass'].fillna("None")再来看一下还有没有缺失值，发现已经没有了，到这我们就完成了对缺失值的处理
all_data_na = (all_data.isnull().sum() / len(all_data)) * 100 all_data_na = all_data_na.drop(all_data_na[all_data_na == 0].index).sort_values(ascending=False) missing_data = pd.DataFrame({'Missing Ratio' :all_data_na}) missing_data.head()此外，根据回归的正态性假设，我们需要检验因变量是否符合正态分布。这里我们对房价的分布进行检测，从房价的直方图分布和qq-plot可以看出，房价的分布是右偏的，需要对其做一些转换让它符合正太分布。
sns.distplot(train['SalePrice'] , fit=norm); (mu, sigma) = norm.fit(train['SalePrice']) print( '\n mu = {:.2f} and sigma = {:.2f}\n'.format(mu, sigma)) plt.legend(['Normal dist. ($\mu=$ {:.2f} and $\sigma=$ {:.2f} )'.format(mu, sigma)], loc='best') plt.ylabel('Frequency') plt.title('SalePrice distribution') fig = plt.figure() res = stats.probplot(train['SalePrice'], plot=plt) plt.show()
我们采用log对数变换对房价进行处理，通过下面转换后房价的分布可以看出，房价转换后符合正态分布。
train["SalePrice"] = np.log1p(train["SalePrice"]) sns.distplot(train['SalePrice'] , fit=norm); (mu, sigma) = norm.fit(train['SalePrice']) print( '\n mu = {:.2f} and sigma = {:.2f}\n'.format(mu, sigma)) plt.legend(['Normal dist. ($\mu=$ {:.2f} and $\sigma=$ {:.2f} )'.format(mu, sigma)], loc='best') plt.ylabel('Frequency') plt.title('SalePrice distribution') fig = plt.figure() res = stats.probplot(train['SalePrice'], plot=plt) plt.show()
对数据做完基本的数据处理后，下面就进入最重要的特征工程，我们需要经过特征预处理、特征抽取和特征选择这三个过程。
特征抽取是通过计算一些统计值或对不同特征进行一些加减乘除变换，获得一些新的特征，增加对模型的估计准确性。我们首先将地下室面积、1楼面积、2楼面积加起来可以得到房屋总面积特征，然后之前可视化时候注意到，建造时间比较近的房子房价比较高，所以新创造一个01特征，如果房屋建造时间在1990年后，则为1，否则是0。
all_data['TotalSF'] = all_data['TotalBsmtSF'] + all_data['1stFlrSF'] + all_data['2ndFlrSF'] all_data['YearBuilt_cut'] = all_data['YearBuilt'].apply(lambda x:1 if x>1990 else 0)特征预处理是将字符串变量进行转化编码，常见的方法有Label Encoding和One_Hot Encoding（独热编码）。首先对于有序性离散变量，我们使用label encoder进行编码
from sklearn.preprocessing import LabelEncoder cols = ('FireplaceQu', 'BsmtQual', 'BsmtCond', 'GarageQual', 'GarageCond', 'ExterQual', 'ExterCond','HeatingQC', 'KitchenQual', 'BsmtFinType1', 'BsmtFinType2', 'Functional', 'BsmtExposure', 'GarageFinish', 'LandSlope', 'LotShape', 'PavedDrive', 'Street', 'CentralAir', 'MSSubClass', 'OverallCond', 'YrSold', 'MoSold') # process columns, apply LabelEncoder to categorical features for c in cols: lbl = LabelEncoder() lbl.fit(list(all_data[c].values)) all_data[c] = lbl.transform(list(all_data[c].values)) # shape print('Shape all_data: {}'.format(all_data.shape))再来看一下数据中的格式，发现之前的英文字符串已经变成了数字
all_data.head(5)据集中还有部分非有序性离散变量，我们通过独热编码将他们转换成哑变量的形式。
all_data = pd.get_dummies(all_data)下面我们进行特征工程的最后一步，特征选择，为避免多重共线性问题，下面我们会识别皮尔森相关性系数大于0.9的特征，并将这些特征删除。
threshold = 0.9 corr_matrix = all_data.corr().abs() corr_matrix.head() upper = corr_matrix.where(np.triu(np.ones(corr_matrix.shape), k=1).astype(np.bool)) upper.head() to_drop = [column for column in upper.columns if any(upper[column] > threshold)] print('There are %d columns to remove.' % (len(to_drop))) all_data = all_data.drop(columns = to_drop) all_data.shape我们将6列相关性系数大于0.9的列删除来避免多重共线性。
对于这次比赛，我们将尝试使用正则化线性回归模型，L1正则 lasso 和L2正则 Ridge。 首先导入一些所需要的包
from sklearn.linear_model import Ridge, RidgeCV, ElasticNet, LassoCV, LassoLarsCV from sklearn.model_selection import cross_val_score然后建立Sklearn的cross-valu-score函数去作为模型准确性的评估
def rmse_cv(model): rmse= np.sqrt(-cross_val_score(model, train, y_train, scoring="neg_mean_squared_error", cv = 5)) return(rmse)接着导入ridge模型，并对模型进行调参，主要需要调节的参数是alpha，alpha决定模型的正则化程度。正则化程度越高，我们的模型越不容易过度拟合。我们列出不同的alpha取值，用我们上面定义的rmse_cv模型评分参数去比较每个参数下模型的分数
model_ridge = Ridge() alphas = [0.05, 0.1, 0.3, 1, 3, 5, 10, 15, 30, 50, 75] cv_ridge = [rmse_cv(Ridge(alpha = alpha)).mean() for alpha in alphas] cv_ridge = pd.Series(cv_ridge, index = alphas) cv_ridge.plot(title = "Validation - Just Do It") plt.xlabel("alpha") plt.ylabel("rmse") cv_ridge
alpha参数选择误差最小的5,然后用训练集对模型进行训练，最后对测试集进行预测，就完成了本项目。
clf = Ridge(alpha=5) clf.fit(train,y_train) predict = clf.predict(test) sub = pd.DataFrame() sub['Id'] = test_ID sub['SalePrice'] = predict sub.head(10) sub.to_csv('submission.csv',index=False)我们将结果上传到Kaggle，得分0.12013，排名前12.18%（565/4637）。
