@程序员，如何在买房时不被宰？
作者 | 胡萝卜
责编 | 胡巍巍
身为程序员，如何自制一个二手房估价模型，以最实惠的价格购得房子？本篇文章讲的就是这件事！不过本文侧重于完整的实现过程和思路，而不是代码部分。
完成后会用自己训练的模型来实战预测下Q房网的成交数据，并且对比下房产大数据平台的“房价网”的估价器。
所用到的Python库：sickit-learn、Numpy、pandas、Matplotlib、beautifulsoup。
机器学习到底是干什么的？
假设你有一位邻居老张，东北人，平时经常去市场买菜，其中买大白菜最多。
由于老张这个人特别会过日子，不希望被摊位随便要价，他就准备自己研究下白菜的价格，买菜都带一个账本记账。买了5次以后记录下来白菜的价格是这样的：
根据账本他推算了下白菜的价格，大约是5元一斤。如果买4斤白菜的话，正常价格应该是20元左右，这个特别简单。
数学上描述这个关系是用一个线性方程y = ax，其中x代表白菜的斤数，y代表白菜的价格。a在这个案例里面等于5（每斤白菜的单价），而y = 5x就是白菜价格的模型。
现在老张觉得这个估价方法挺管用，由于价格猜测得准确也没有摊位敢胡乱开价。太平过了一段时间后，有一天他又去买菜，这次摊主却告诉他4斤白菜要28块了。老张很疑惑，问摊主原因。
摊主说，哦，因为这批白菜比较新鲜，刚采摘下来不到24小时，所以卖的比较贵。而且，由于家里要买一块新菜地，为了增加点收入，以后白菜也都要按照新鲜程度来卖了…...
这样一来，原来那个模型就没用了。于是老张又记录了一段时间的账本，与之前不同的是他这次还记录了白菜的新鲜程度。这个账本是这样的：
这么一看还是能隐约知道白菜的价格与斤数有关，但是似乎关系不那么明显了。现在老张怎么来估算白菜的价格？
其实白菜的价格分布还是有规律的。假设把老张的账本投射到一个三维空间里是这样的：
可以看到几乎所有的点在三维空间中都处于同一个平面范围上。根据高中数学我们会知道三维空间平面的表达式是z = ax + by + c。
这里的a不再代表白菜的单价了，仅仅是白菜重量的系数。新的账本的表达式是z = 7x - 3.5y + 3.5。
这个白菜价格估算方法用的就是机器学习中最基础也最经典的一种算法，叫做Linear Regression线性回归。
老张的账本在机器学习中叫做训练数据集，重量、新鲜程度等描述白菜属性的数值或者分类叫做特征。
而求解模型中a,b,c参数具体数值，使得它对所有预测结果与真实值之间综合误差最小的过程就叫做模型的拟合。
这看起来很OK，可是你一定会说现实生活中的问题根本不是这样的啊！事实上老张可能会去多个摊位买菜，每次去哪个摊位都是随机的。
又或者白菜有大有小，白菜大小也会影响价格。白菜也可能在冬天更便宜，夏天更贵等等…...
说的没错，实际生活中我们往往需要多个特征来描述问题，对应一个多维以上的空间。
数据分布也并非总是线性的，可能是一个曲面或者高维超曲面，数据也可能并不会正好都在某个曲面上等等。
这样问题就来了，大多数的人对高维想象都很难，怎么去解一个高维空间的问题哪？
这个时候计算机就发挥作用了。因为在机器“眼”里世界是数学抽象的，它不需要理解或者想象高维空间，只需要将低维空间的运算规则推广到高维空间即能处理一系列求解。
进一步，计算机寻找最佳参数的方式叫做梯度下降法（理解需要一点入门级别的微积分。)
这种方法很容易让我联想起曾经看过的一系列与时间循环有关的科幻电影。
电影讲述主角被困于某段时间线中，每次死亡以后又都会回到时间线的最初。
每次重生到再死亡的过程中也都会多获得关于事件一点的信息，最后所有拼凑起来的信息还原了事件的真相。
如果把电影中的每次重生到死亡的过程看成一次“迭代”，那么计算机寻找最佳参数的梯度下降法就是迭代，每次迭代都向最优方向前进一点点，当迭代非常多的次数后最终就能非常逼近最优参数。
这跟现实生活中人的学习方法很接近，所以机器学习叫Machine Learning而不叫Machine Fitting，或者Machine Predicting什么的…...
区别只是计算机的迭代时间可能1万次只用了几秒钟，而人类，因为现实中不会真的有时间循环让你去重复经历同一件事，可能在重大问题上迭代个几次大概一辈子就过去了…...
所以，概括地说机器学习做的事情就是输入训练数据集，给定一种建模方式，计算机自动寻找最佳拟合参数使模型可以描述数据集中输入和输出的对应关系。并用这个模型来预测新输入数据的过程。
准备工作
既然已经知道机器学习是什么了，我们就要着手开始制作自己的模型了。参考上方关系图，我们需要准备点什么哪？
首先，我们需要一个开源库，不用自己写一大堆晦涩艰深的数学公式去指导机器计算，只要传参数就可以傻瓜式操作了。
Python有一个机器学习的库Scikit-learn就很好用。为了熟悉库，需要看下使用文档。
我们还需要一台能计算的电脑，我就用了公司配发的低配行政笔记本电脑 。
解决特定领域问题的时候，该领域的专业知识会帮助你。比如，你要通过人脸表情照片去识别笑容，就需要了解一点图形学，知道计算机“看”照片是一个像素矩阵，每个像素点的灰度值是一个数字等等。
最后，就是最重要的数据源了。选择数据的质量和规模是直接影响模型表现的最重要因素。更多的时候可能我们想的到要解决哪些问题，却根本不知道从哪儿去找数据源…...
这里我们选择做一个上海二手房的成交估价模型，因为相对数据更好采集。
采集还是通过爬虫来实现，对象则是最受广大爬虫玩家欢迎的房产网站“链家网”。
动手采集前我们需要先看下链家“二手房成交”板块房产详情页，分析下大致哪些特征可能对判断成交价格有用。
区域，板块、小区名称、成交价格、成交日期这几项是必须采集的。挂牌价格、成交天数、带看、关注、浏览量这几项假设想进一步分析成交时间的话会有用，可采可不采。
户型、楼层、面积、朝向、梯户比这几项是直觉与价格有关的因素，所以采下来。
小区信息中建筑年份、物业费、总楼栋数和总户数这四个特征我们也认为与成交价格有关，所以采集下来。
这里你可能会问为啥没有采集“小区均价”哪？估算房价最直接的不应该是小区均价吗？
其实是因为链家网站推算小区均价的逻辑，这里的小区均价计算的是“挂牌价格”的均价。
留意下案例中这套房子，成交均价是46259元/平，并且成交时间就是离现在很近的9.30日。
而小区的平均挂牌均价是57507元/平。直觉告诉我们房产虽然具有投资属性但并不可能在20天内有这样大规模的波动，既然我们研究的问题是成交，那么就以成交价格为准。
最后要采集的就是配套了。作为一个实验案例我这里并没有采集医院，学校等信息。而是着重采集了小区经纬度和周边1.5公里直线距离内的地铁站个数，地铁线路条数。
整个爬虫的代码是比较简单的，类似爬取“链家网”的博文CSDN上可以找到很多，用到的库就是beautifulsoup，这边就不赘述了。
贴一下调试爬取地铁配套部分的代码吧，这里需要调用下百度地图的API来定位到小区经纬度，并且用POI来查找周边地铁站个数和地铁线路数，返回json格式再解析出来。
# 输入上海任意小区名，打印出周边1.5公里直接距离内的地铁站个数，名称，地铁线路和步行距离
name_estate = input("输入小区名字： ")
#中文转码utf-8
name_estate_quoted = quote("上海"+name_estate)
#调用百度地图API获得小区经纬度
find_location = "http://api.map.baidu.com/geocoder/v2/?address="+name_estate_quoted+"&output=json&ak=你申请的ak"
page = urlopen(find_location).read()
content = json.loads(page,encoding= "utf-8")
o_lat = content["result"]["location"]["lat"]
o_lng = content["result"]["location"]["lng"]
#用获取的小区经纬度作为参数，调用百度地图的POI查找周边地铁站。radius传范围大小（米），output返回格式
find_metro = r"http://api.map.baidu.com/place/v2/search?query=%E5%9C%B0%E9%93%81&location="+str(o_lat)+","+str(o_lng)+"&radius=1500&output=json&ak=你申请的ak"
page = urlopen(find_metro).read()
results = json.loads(page,encoding = "utf-8")["results"]
counts = len(results)
print("【"+name_estate+"】周边1.5公里范围内共有"+str(counts)+"个地铁站，分别是：")
for result in results:
d_lat = result["location"]["lat"]
d_lng = result["location"]["lng"]
#调取百度地图算路,mode传出行方式（步行，驾车等), origin传出发地经纬度，destination传目的地经纬度
calcu = "http://api.map.baidu.com/direction/v1?mode=walking&origin="+str(o_lat)+","+str(o_lng)+"&destination="+str(d_lat)+","+str(d_lng)+r"&origin_region=%E4%B8%8A%E6%B5%B7&destination_region=%E4%B8%8A%E6%B5%B7&output=json&ak=你申请的ak"
obj = urlopen(calcu).read()
content = json.loads(obj,encoding = "utf-8")
res = content["result"]["routes"][0]
distance = res["distance"]
duration = res["duration"]
print(result["name"]+" "+result["address"]+" 步行距离约"+str(int(distance))+"米"+" 耗时约"+str(int(duration/60))+"分钟")
另外一个采集前需要考虑的问题是，我们是否有必要控制下数据时效性？
假设我们打算只估计近期二手房成交价格，那么因为价格的波动，太久远的数据反而可能让模型产生偏差。
所以我们圈定了一个时间范围为7月至今。最后，采集完成后就得到了大约7901组数据。
对数据的清洗和预处理
到这里为止上面那堆数据还不能直接拿来训练模型，我们还需要对其进行清洗和预处理。
▌处理虚拟变量
第一个问题是机器无法处理像类似“两梯三户”这种文字特征，或者说这种表述方式无法给予机器有效信息。
一种处理方法是我们将这个特征做成“虚拟变量”或这叫One-Hot编码，其实就是一个01矩阵。
打个比方来说，在梯户比这个特征上假设可能出现的结果有“一梯两户”，“一梯四户”，“两梯三户”，“两梯四户”这4种可能性，一个“一梯四户”房产就表示为下面这种形式
这个编码可以用pandas的get_dummies方法来实现，非常方便。假设你不想逐一设定列名的话，使用get_dummies之前唯一要小心的点在于要确认所有数值型的数据类型不是object类型，否则get_dummies是会把数值类型特征也虚拟变量化的。
因为虚拟变量会大大增大特征维度，造成计算量上升。而梯户比的实际含义是数值，也可以直接处理成两列，一列代表梯数，一列代表户数。显然“梯户比”这个特征这里处理成数值更好。
最终我们直接去除了“小区信息”，没有把它作为输入变量。原因一是假设对小区进行虚拟变量变换的话会大大增加数据维度从而对计算性能提出更高的要求。
二是我们目前的数据量没有足够大到覆盖上海所有小区，假设预测新数据的小区并没有出现训练数据集里则会造成特征不一致的问题，代码会直接报错。
▌填充缺失值
其次，网络采集的数据都可能会存在大量缺失值。比如下面这种“暂无信息”。
在我们这个训练数据集里，有缺失值的数据有703条，几乎占了总数据量的9%。如果我们不想损失掉这些数据就不能粗暴的将它们删除，而是要设定一定的方式对缺失值填空。
这里我们可以用numpy的isnull方法来查找下哪些列有缺失值，发现是“成交时间”、“朝向”、“电梯”、“建设年份”和“物业费”这5列。
#处理“暂无信息”或者“暂无数据”
for_training = for_training.apply(lambda x:x.replace("暂无信息",np.nan).replace("暂无数据",np.nan)
for_training.isnull().any()
其中成交天数我们最终不打算把它作为输入特征，可以随便给它一个值后不用管它。朝向我们统一给它填充“南”，有无电梯我们按照2000年前<=6层的建筑估算“无”，其余估算“有”来填充。
建设年份按照同板块楼盘建设年份的平均数来估。物业费则按照同建设年份物业费的平均数来估。
▌查找异常值
上面这些都完成了以后还需要观察下现有数据。
想象下在老张买菜的案例里面，如果他记录账本的那段时间正好碰到白菜大减价，那么输入大量减价后的价格特征，模型一定会产生偏斜。
在二手房的问题上像下面这种成交价格低的不可思议的（相对上海房价来说），或者挂牌价格和成交价格相差巨大的，就可以判定为典型异常值。
这里我们用统计学的分箱图来排除异常值，我们计算下成交均价的log变换后做下分箱：
分箱图的看法是这样的，中间红线代表“中位数”，箱体的上下边缘分别是“上四分位”和“下四分位”。上下四分位间的距离叫做“四分位距”。而上下超过1.5倍四分位距的数值都被判断为异常值。这里大约要删除53组数据。
删除后可以看到二手房成交均价的分布1.598~12.612万之间，较为符合我们对上海房价的逻辑常识认知了。
完成这步后，最后得到了一个7833x245的数据集。去除不作为输入的信息，基本上可以知道我们输入数据的维度在240左右。
训练模型、调参和可视化
我们来为模型选择一种算法，这里预测二手房成交价格是个回归问题，我们选择RandomForestRegression随机森林回归。
与一开始老张买菜的案例不同，二手房问题的复杂度高的多。线性模型我在这里也调试了下，表现最好的情况是L1正则化以后的Lasso可以达到0.84分（满分为1，表示100%的数据可用模型解释），这个分数不算太低。
但树集成类算法在这个问题上可以表现更好。关于随机森林的原理有兴趣的可以自行百度，简单来说可以理解为N棵随机的决策树通过分叉后覆盖所有数据，然后再取平均。
因为scikit-learn是个傻瓜式工具包，我们只需要为算法调节一些参数。分别是随机树的棵树（n_estimators）和树的最大深度（max_depth）。在scikit-learn里面最佳参数的查找也是可以用网格搜索grid_search查找的。
#读取清洗好的数据集
data = pd.read_csv(r"你的目录\shhouse_dummies.csv",header = 0,encoding = "gbk")
#打乱数据集
data = data.reindex(np.random.permutation(data.index))
#设计成交价格为预测目标
target = data["deal_price"]
#删除不作为输入特征的列
data.drop("deal_price",axis = 1,inplace = True)
data.drop("post_price",axis = 1,inplace = True)
data.drop("deal_days",axis = 1,inplace = True)
data.drop("price_per_area",axis = 1,inplace = True)
data.drop("community",axis = 1,inplace = True)
#分割数据（注：正规做法是这里是要将数据集分割为训练集和测试集的，由于我们下面会启动五折交叉验证，为了节省数据集就不再分割了）
#X_train,X_test,y_train,y_test = train_test_split(data,target,random_state = 1)
X_train = data
y_train = target
#调用scikit-learn的网格搜索，传入参数选择范围，并且制定随机森林回归算法，cv = 5表示5折交叉验证
param_grid = {"n_estimators":[5,10,50,100,200,500],"max_depth":[5,10,50,100,200,500]}
grid_search = GridSearchCV(RandomForestRegressor(),param_grid,cv = 5)
#让模型对训练集和结果进行拟合
grid_search.fit(X_train,y_train)
print(np.around(grid_search.best_score_,2))
我们尝试两个参数在[5,10,50,100,200,500]中各种排列组合的可能性，并对训练集进行5折交叉验证（平均分成五分，每次各用不同的四份来训练，用剩下的一份来测试）来选出最优参数。
完了以后运行代码就是等待了。根据机器的计算性能需要等待不同的时间，我的行政笔本等待的时间约为20-30分钟左右。
结束后可以看到最终我们获得了一个约0.90分的模型，即约90%的数据可以用模型来解释，这高于了线性模型约6个百分点。该模型最佳的参数选择是500棵树，50层深度。
我们还可以将不同参数的组合结果用Matplotlib的imshow可视化一下，代码如下：
fig = plt.figure(figsize = (16,9))
ax = fig.add_subplot(1,1,1,facecolor = "whitesmoke",alpha = 0.2)
ax.imshow(df,cmap = "summer")
ax.set_xlim(-0.5,5.5)
ax.set_ylim(-0.5,5.5)
ax.set_xticklabels([0,5,10,50,100,200,500],fontsize = 18)
ax.set_yticklabels([0,5,10,50,100,200,500],fontsize = 18)
ax.set_xlabel("n_estimators",fontsize = 18)
ax.set_ylabel("max_depth",fontsize = 18)
for i in range(0,6,1):
for j in range(0,6,1):
ax.text(i,j,str(np.around(df.iloc[j,i],3)),fontsize = 15,verticalalignment="center",horizontalalignment="center",color = "black")
▌得到如下结果：
这就完了？你可能会说这样一点都不直观啊！我该怎么去解释这个完成的模型是什么样的哪？
还好，scikit-learn的树算法还提供了一个叫特征权重的属性。我们可以把这个属性调出来可视化一下，看下从机器的“眼睛”如何解读影响房价的这些特征因素。代码是这样的：
#特征重要前十位性可视化
features = X_train.columns
importance = grid_search.best_estimator_.feature_importances_
fi = pd.Series(importance,index = features)
fi = fi.sort_values(ascending = False)
ten = fi[:10]
fig = plt.figure(figsize = (16,9))
ax = fig.add_subplot(1,1,1,facecolor = "whitesmoke",alpha = 0.2)
ax.grid(color = "grey",linestyle=":",alpha = 0.8,axis = "y")
ax.barh(ten.index,ten.values,color = "dodgerblue")
ax.set_xticklabels([0.0,0.1,0.2,0.3,0.4,0.5,0.6],fontsize = 22)
ax.set_yticklabels(ten.index,fontsize = 22)
ax.set_xlabel("importance",fontsize = 22)
▌结果如下：
上图展示了机器评价重要程度前10位的特征，所有重要程度的和为1。
可以看出排在第1位的是“面积”，的确符合常识，面积是与总价关联性最强的因素，影响权重在0.6左右。
第2位有点出人意料，机器选中的是“物业费”。一种可能性是物业费的高低反应了小区的档次，小区的档次是影响不同楼盘单价差异的重要因素。
第3,4位分别是“地铁站数量”和“地铁线路数”，这个比较符合常识，交通方便理论上房价就会提升。
第5,6位是经纬度，也就是小区在上海的实际位置，反应的是区位。
后7,8,9,10位都是与小区有关的特征，有“总户数（单元数）”、“建设年份”、“总层数（高度）”和“总户数”。
排除建设年份这点，楼盘设计师们应该会很高兴，因为机器认为一个楼盘的规划设计参数的确影响了房价。
模型预测实战VS“房价网”估价器
坚持看到这里的你一定希望看下模型的应用效果吧？这肯定比给出一个0.90分的评分更直观，我们就来试一下。
由于我们的训练数据集来自链家，测试的时候就不能再用链家数据来测试，某则模型会给你一个100%准确的预测结果。我们需要一个全新的数据集。
这里我们选择了也可以查到成交数据的“Q房网”。选取Q房第一页8月底至今的成交数据20条。
为了让测试更有意思一点，我们特别让模型对比了下“房产大数据平台”的“房价网”的估价器。
下面就是最终结果。
自制模型大多数情况误差都小于12%，总体要好于“房价网”的估价器，其中高于15%的误差总共出现了3次，在预测区域高均价房产上的表现较差。而“房价网”估价器高于15%的误差出现了6次，总体误差范围更大。
查找误差原因和改进模型
为了进一步改进，尝试找下自制模型较大误差部分产生的原因。
以第7号数据“大华铂金华府”，误差27.26%为例，我们怀疑可能原因是数据分布造成的。为了验证用Matplotlib来看下该楼盘所在的宝山区-大场板块训练数据集的成交均价分布。
整个宝山区的价格分布集中于3-5万之间，“大华铂金华府”在整个宝山区属于偏右侧尾部均价偏高的楼盘，成交均价约为62222元/平方米。
而它所在大场板块，训练数据集中最高成交均价也仅为4.897万。果然与我们猜测的情况是一致的。这里也可以看出我们开头提到的模型的性能表现主要取决于训练集的数据规模和质量。
优化这个模型的一种方式（可能）是放宽时间维度以换取更大的训练数据集体量，或者多渠道获取数据集，以保证训练集含有一定数量的高成交均价样本。
作者：胡萝卜，CSDN博客专家。普通的地产从业人员，职业是市场BD。自学Python和数据分析，希望让技术接地气，解决实际工作和日常生活中的问题。
微信改版了，
想快速看到CSDN的热乎文章，
赶快把CSDN公众号设为星标吧，
打开公众号，点击“设为星标”就可以啦！
“征稿啦”
CSDN 公众号秉持着「与千万技术人共成长」理念，不仅以「极客头条」、「畅言」栏目在第一时间以技术人的独特视角描述技术人关心的行业焦点事件，更有「技术头条」专栏，深度解读行业内的热门技术与场景应用，让所有的开发者紧跟技术潮流，保持警醒的技术嗅觉，对行业趋势、技术有更为全面的认知。
如果你有优质的文章，或是行业热点事件、技术趋势的真知灼见，或是深度的应用实践、场景方案等的新见解，欢迎联系 CSDN 投稿，联系方式：微信（guorui_1118，请备注投稿+姓名+公司职位），邮箱（guorui@csdn.net）。
推荐阅读：