(干貨)數據分析案例--以上海二手房為例

如果你手上有一批數據，你可能應用統計學、挖掘算法、可視化方法等技術玩轉你的數據，但你沒有數據的時候，我該怎么玩呢？接下來就帶著大家玩玩沒有數據情況下的數據分析。
本文從如下幾個目錄詳細講解數據分析的流程：
1、數據源的獲??；
2、數據探索與清洗；
3、模型構建（聚類算法和線性回歸）；
4、模型預測；
5、模型評估；
一、數據源的獲取
正如本文的題目一樣，我要分析的是上海二手房數據，我想看看哪些因素會影響房價？哪些房源可以歸為一類？我該如何預測二手房的價格？可我手上沒有這樣的數據樣本，我該如何回答上面的問題呢？
互聯網時代，網絡信息那么發達，信息量那么龐大，隨便找點數據就夠喝一壺了。前幾期我們已經講過了如何從互聯網中抓取信息，采用Python這個靈活而便捷的工具完成爬蟲，例如：
通過Python抓取天貓評論數據
使用Python實現豆瓣閱讀書籍信息的獲取
使用Python爬取網頁圖片
當然，上海二手房的數據仍然是通過爬蟲獲取的，爬取的平臺來自于鏈家，頁面是這樣的：

我所需要抓取下來的數據就是紅框中的內容，即上海各個區域下每套二手房的 小區名稱、戶型、面積、所屬區域、樓層、朝向、售價及單價 。先截幾張Python爬蟲的代碼，源代碼和數據分析代碼寫在文后的鏈接中，如需下載可以到指定的百度云盤鏈接中下載。

上面圖中的代碼是構造所有需要爬蟲的鏈接。

上面圖中的代碼是爬取指定字段的內容。
爬下來的數據是長這樣的（總共28000多套二手房）：

二、數據探索與清洗（一下均以R語言實現）
當數據抓下來后，按照慣例，需要對數據做一個探索性分析，即了解我的數據都長成什么樣子。
1、戶型分布
# 戶型分布
library(ggplot2)
type_freq <- data.frame(table(house$戶型))
# 繪圖
type_p <- ggplot(data = type_freq, mapping = aes(x = reorder(Var1, -Freq),y = Freq)) + geom_bar(stat = 'identity', fill = 'steelblue') + theme(axis.text.x  = element_text(angle = 30, vjust = 0.5)) + xlab('戶型') + ylab('套數')
type_p

我們發現只有少數幾種的戶型數量比較多，其余的都非常少，明顯 屬于長尾分布類型（嚴重偏態） ，所以，考慮將1000套一下的戶型統統歸為一類。
# 把低于一千套的房型設置為其他
type <- c('2室2廳','2室1廳','3室2廳','1室1廳','3室1廳','4室2廳','1室0廳','2室0廳')
house$type.new <- ifelse(house$戶型 %in% type, house$戶型,'其他')
type_freq <- data.frame(table(house$type.new))
# 繪圖
type_p <- ggplot(data = type_freq, mapping = aes(x = reorder(Var1, -Freq),y = Freq)) + geom_bar(stat = 'identity', fill = 'steelblue') + theme(axis.text.x  = element_text(angle = 30, vjust = 0.5)) + xlab('戶型') + ylab('套數')
type_p

2、二手房的面積和房價的分布
# 面積的正態性檢驗
norm.test(house$面積)

# 房價的正態性檢驗
norm.test(house$價格.W.)

上面的norm.test函數是我自定義的函數 ，函數代碼也在下文的鏈接中，可自行下載。從上圖可知， 二手房的面積和價格均不滿足正態分布，那么就不能直接對這樣的數據進行方差分析或構建線性回歸模型 ，因為這兩種統計方法，都要求正態性分布的前提假設，后面我們會將講解如何處理這樣的問題。
3、二手房的樓層分布
原始數據中關于樓層這一變量，總共有151種水平，如地上5層、低區/6層、中區/11層、高區/40層等，我們覺得有必要將這151種水平設置為低區、中區和高區三種水平，這樣做有助于后面建模的需要。
# 把樓層分為低區、中區和高區三種
house$floow <- ifelse(substring(house$樓層,1,2) %in% c('低區','中區','高區'), substring(house$樓層,1,2),'低區')
# 各樓層類型百分比分布
percent <- paste(round(prop.table(table(house$floow))*100,2),'%',sep = '')
df <- data.frame(table(house$floow))
df <- cbind(df, percent)
df

可見，三種樓層的分布大體相當，最多的為高區，占了36.1%。
4、上海各地區二手房的均價
# 上海各區房價均價
avg_price <- aggregate(house$單價.平方米., by = list(house$區域), mean)
#繪圖
p <- ggplot(data = avg_price, mapping = aes(x = reorder(Group.1, -x), y = x, group = 1)) + geom_area(fill = 'lightgreen') + geom_line(colour = 'steelblue', size = 2) + geom_point() + xlab('') + ylab('均價')
p

很明顯，上海二手房價格最高的三個地區為：靜安、黃埔和徐匯，均價都在7.5W以上，價格最低的三個地區為：崇明、金山和奉賢。
5、房屋建筑時間缺失嚴重

建筑時間這個變量有6216個 缺失，占了總樣本的22% 。雖然缺失嚴重，但我也不能簡單粗暴的把該變量扔掉， 所以考慮到按各個區域分組，實現眾數替補法 。這里構建了兩個自定義函數：
library(Hmisc)
# 自定義眾數函數
stat.mode <- function(x, rm.na = TRUE){
  if (rm.na == TRUE){
    y = x[!is.na(x)]
  }
  res = names(table(y))[which.max(table(y))]
  return(res)
}
# 自定義函數，實現分組替補
my.impute <- function(data, category.col = NULL,
                      miss.col = NULL, method = stat.mode){
  impute.data = NULL
  for(i in as.character(unique(data[,category.col]))){
    sub.data = subset(data, data[,category.col] == i)
    sub.data[,miss.col] = impute(sub.data[,miss.col], method)
    impute.data = c(impute.data, sub.data[,miss.col])
  }
  data[,miss.col] = impute.data
  return(data)
}

# 將建筑時間中空白字符串轉換為缺失值
house$建筑時間[house$建筑時間 == ''] <- NA
#分組替補缺失值，并對數據集進行變量篩選

final_house <- subset(my.impute(house, '區域', '建筑時間'),select = c(type.new,floow,面積,價格.W.,單價.平方米.,建筑時間))
#構建新字段，即建筑時間與當前2016年的時長
final_house <- transform(final_house, builtdate2now = 2016-as.integer(substring(as.character(建筑時間),1,4)))
#刪除原始的建筑時間這一字段
final_house <- subset(final_house, select = -建筑時間)
最終完成的干凈數據集如下：

接下來就可以針對這樣的干凈數據集，作進一步的分析，如聚類、線性回歸等。
三、模型構建
這么多的房子，我該如何把它們分分類呢？即應該把哪些房源歸為一類？這就要用到聚類算法了， 我們就使用簡單而快捷的k-means算法實現聚類的工作 。但聚類前，我需要掂量一下我該聚為幾類？根據 聚類原則：組內差距要小，組間差距要大 。我們繪制不同類簇下的組內離差平方和圖，聚類過程中，我們選擇面積、房價和單價三個數值型變量：
tot.wssplot <- function(data, nc, seed=1234){
  #假設分為一組時的總的離差平方和              
  tot.wss <- (nrow(data)-1)*sum(apply(data,2,var))
  for (i in 2:nc){
    #必須指定隨機種子數
    set.seed(seed)
    tot.wss[i] <- kmeans(data, centers=i, iter.max = 100)$tot.withinss
  }
  plot(1:nc, tot.wss, type="b", xlab="Number of Clusters",
       ylab="Within groups sum of squares",col = 'blue',
       lwd = 2, main = 'Choose best Clusters')
}
# 繪制不同聚類數目下的組內離差平方和
standrad <- data.frame(scale(final_house[,c('面積','價格.W.','單價.平方米.')]))
myplot <- tot.wssplot(standrad, nc = 15)

當把所有樣本當作一類時，離差平方和達到最大，隨著聚類數量的增加，組內離差平方和會逐漸降低，直到極端情況，每一個樣本作為一類，此時組內離差平方和為0。從上圖看，聚類數量在5次以上，組內離差平方降低非常緩慢，可以把拐點當作5，即聚為5類。
# 將樣本數據聚為5類
set.seed(1234)
clust <- kmeans(x = standrad, centers = 5, iter.max = 100)
table(clust$cluster)

# 按照聚類的結果，查看各類中的區域分布
table(final_house$區域,clust$cluster)

# 各戶型的平均面積
aggregate(final_house$面積, list(final_house$type.new), mean)

# 按聚類結果，比較各類中房子的平均面積、平均價格和平均單價
aggregate(final_house[,3:5], list(clust$cluster), mean)

從平均水平來看，我大體可以將28000多套房源合成為如下幾種說法：
a、大戶型（3室2廳、4室2廳），屬于第2類。 平均面積都在130平以上，這種大戶型的房源主要分布在青浦、黃埔、松江等地（具體可從 各類中的區域分布圖可知 ）。
b、地段型（房價高），屬于第1類。 典型的區域有黃埔、徐匯、長寧、浦東等地（具體可從 各類中的區域分布圖可知 ）。
c、大眾蝸居型（面積小、價格適中、房源多），屬于第4和5類。 典型的區域有寶山、虹口、閔行、浦東、普陀、楊浦等地
d、徘徊型（大戶型與地段型之間的房源），屬于第3類。 典型的區域有奉賢、嘉定、青浦、松江等地。這些地區也是將來迅速崛起的地方。
# 繪制面積與單價的散點圖，并按聚類進行劃分
p <- ggplot(data = final_house[,3:5], mapping = aes(x = 面積,y = 單價.平方米., color = factor(clust$cluster)))
p <- p + geom_point(pch = 20, size = 3)
p + scale_colour_manual(values = c("red","blue", "green", "black", "orange"))

接下來我想 借助于已有的數據（房價、面積、單價、樓層、戶型、建筑時長、聚類水平）構建線性回歸方程 ，用于房價因素的判斷及預測。由于數據中有離散變量，如戶型、樓層等，這些變量入模的話需要對其 進行啞變量處理 。
# 構造樓層和聚類結果的啞變量
# 將幾個離散變量轉換為因子，目的便于下面一次性處理啞變量
final_house$cluster <- factor(clust$cluster)
final_house$floow <- factor(final_house$floow)
final_house$type.new <- factor(final_house$type.new)
# 篩選出所有因子型變量
factors <- names(final_house)[sapply(final_house, class) == 'factor']
# 將因子型變量轉換成公式formula的右半邊形式
formula <- f <- as.formula(paste('~', paste(factors, collapse = '+')))
dummy <- dummyVars(formula = formula, data = final_house)
pred <- predict(dummy, newdata = final_house)
head(pred)

# 將啞變量規整到final_house數據集中
final_house2 <- cbind(final_house,pred)
# 篩選出需要建模的數據
model.data <- subset(final_house2,select = -c(1,2,3,8,17,18,24))
# 直接對數據進行線性回歸建模
fit1 <- lm(價格.W. ~ .,data = model.data)
summary(fit1)

從體看上去還行，只有建筑時長和2室0廳的房型參數不顯著，其他均在0.01置信水平下顯著。 不要贊贊自喜 ，我們說，使用線性回歸是有假設前提的，即因變量滿足正態或近似于正態分布，前面說過，房價明顯在樣本中是偏態的，并不服從正態分布，所以這里 使用COX-BOX變換處理 。根據COX-BOX變換的lambda結果，我們針對y變量進行轉換，即：
# Cox-Box轉換
library(car)
powerTransform(fit1)
 
根據結果顯示， 0.23非常接近上表中的0值，故考慮將二手房的價格進行對數變換。
fit2 <- lm(log(價格.W.) ~ .,data = model.data)
summary(fit2)

這次的結果就明顯比fit1好很多 ，僅有樓層的中區在0.1置信水平下顯著，其余變量均在0.01置信水平下顯著， 而且調整的R方值也提高到了94.3% ，即這些自變量對房價的解釋度達到了94.3%。
最后我們再看一下，關于最終模型的診斷結果：
# 使用plot方法完成模型定性的診斷
opar <- par(no.readonly = TRUE)
par(mfrow = c(2,2))
plot(fit2)
par(opar)

從上圖看，基本上滿足了線性回歸模型的幾個假設，即：殘差項服從均值為0（左上），標準差為常數（左下）的正態分布分布（右上）?；谶@樣的模型，我們就可以有針對性的預測房價啦~