House_price
完整的一次数据预测项目，使用的是ames房价预测数据集。首先查看数据
str(data)## 'data.frame': 2930 obs. of 83 variables:
## $ PID : int 526301100 526350040 526351010 526353030 527105010 527105030 527127150 527145080 527146030 527162130 ...
## $ MS_SubClass : chr "One_Story_1946_and_Newer_All_Styles" "One_Story_1946_and_Newer_All_Styles" "One_Story_1946_and_Newer_All_Styles" "One_Story_1946_and_Newer_All_Styles" ...
## $ MS_Zoning : chr "Residential_Low_Density" "Residential_High_Density" "Residential_Low_Density" "Residential_Low_Density" ...
## $ Lot_Frontage : int 141 80 81 93 74 78 41 43 39 60 ...
## $ Lot_Area : int 31770 11622 14267 11160 13830 9978 4920 5005 5389 7500 ...
## $ Street : chr "Pave" "Pave" "Pave" "Pave" ...
## $ Alley : chr "No_Alley_Access" "No_Alley_Access" "No_Alley_Access" "No_Alley_Access" ...
## $ Lot_Shape : chr "Slightly_Irregular" "Regular" "Slightly_Irregular" "Regular" ...
## $ Land_Contour : chr "Lvl" "Lvl" "Lvl" "Lvl" ...
## $ Utilities : chr "AllPub" "AllPub" "AllPub" "AllPub" ...
## $ Lot_Config : chr "Corner" "Inside" "Corner" "Corner" ...
## $ Land_Slope : chr "Gtl" "Gtl" "Gtl" "Gtl" ...
## $ Neighborhood : chr "North_Ames" "North_Ames" "North_Ames" "North_Ames" ...
## $ Condition_1 : chr "Norm" "Feedr" "Norm" "Norm" ...
## $ Condition_2 : chr "Norm" "Norm" "Norm" "Norm" ...
## $ Bldg_Type : chr "OneFam" "OneFam" "OneFam" "OneFam" ...
## $ House_Style : chr "One_Story" "One_Story" "One_Story" "One_Story" ...
## $ Overall_Qual : chr "Above_Average" "Average" "Above_Average" "Good" ...
## $ Overall_Cond : chr "Average" "Above_Average" "Above_Average" "Average" ...
## $ Year_Built : int 1960 1961 1958 1968 1997 1998 2001 1992 1995 1999 ...
## $ Year_Remod_Add : int 1960 1961 1958 1968 1998 1998 2001 1992 1996 1999 ...
## $ Roof_Style : chr "Hip" "Gable" "Hip" "Hip" ...
## $ Roof_Matl : chr "CompShg" "CompShg" "CompShg" "CompShg" ...
## $ Exterior_1st : chr "BrkFace" "VinylSd" "Wd Sdng" "BrkFace" ...
## $ Exterior_2nd : chr "Plywood" "VinylSd" "Wd Sdng" "BrkFace" ...
## $ Mas_Vnr_Type : chr "Stone" "None" "BrkFace" "None" ...
## $ Mas_Vnr_Area : int 112 0 108 0 0 20 0 0 0 0 ...
## $ Exter_Qual : chr "Typical" "Typical" "Typical" "Good" ...
## $ Exter_Cond : chr "Typical" "Typical" "Typical" "Typical" ...
## $ Foundation : chr "CBlock" "CBlock" "CBlock" "CBlock" ...
## $ Bsmt_Qual : chr "Typical" "Typical" "Typical" "Typical" ...
## $ Bsmt_Cond : chr "Good" "Typical" "Typical" "Typical" ...
## $ Bsmt_Exposure : chr "Gd" "No" "No" "No" ...
## $ BsmtFin_Type_1 : chr "BLQ" "Rec" "ALQ" "ALQ" ...
## $ BsmtFin_SF_1 : int 2 6 1 1 3 3 3 1 3 7 ...
## $ BsmtFin_Type_2 : chr "Unf" "LwQ" "Unf" "Unf" ...
## $ BsmtFin_SF_2 : int 0 144 0 0 0 0 0 0 0 0 ...
## $ Bsmt_Unf_SF : int 441 270 406 1045 137 324 722 1017 415 994 ...
## $ Total_Bsmt_SF : int 1080 882 1329 2110 928 926 1338 1280 1595 994 ...
## $ Heating : chr "GasA" "GasA" "GasA" "GasA" ...
## $ Heating_QC : chr "Fair" "Typical" "Typical" "Excellent" ...
## $ Central_Air : chr "Y" "Y" "Y" "Y" ...
## $ Electrical : chr "SBrkr" "SBrkr" "SBrkr" "SBrkr" ...
## $ First_Flr_SF : int 1656 896 1329 2110 928 926 1338 1280 1616 1028 ...
## $ Second_Flr_SF : int 0 0 0 0 701 678 0 0 0 776 ...
## $ Low_Qual_Fin_SF : int 0 0 0 0 0 0 0 0 0 0 ...
## $ Gr_Liv_Area : int 1656 896 1329 2110 1629 1604 1338 1280 1616 1804 ...
## $ Bsmt_Full_Bath : int 1 0 0 1 0 0 1 0 1 0 ...
## $ Bsmt_Half_Bath : int 0 0 0 0 0 0 0 0 0 0 ...
## $ Full_Bath : int 1 1 1 2 2 2 2 2 2 2 ...
## $ Half_Bath : int 0 0 1 1 1 1 0 0 0 1 ...
## $ Bedroom_AbvGr : int 3 2 3 3 3 3 2 2 2 3 ...
## $ Kitchen_AbvGr : int 1 1 1 1 1 1 1 1 1 1 ...
## $ Kitchen_Qual : chr "Typical" "Typical" "Good" "Excellent" ...
## $ TotRms_AbvGrd : int 7 5 6 8 6 7 6 5 5 7 ...
## $ Functional : chr "Typ" "Typ" "Typ" "Typ" ...
## $ Fireplaces : int 2 0 0 2 1 1 0 0 1 1 ...
## $ Fireplace_Qu : chr "Good" "No_Fireplace" "No_Fireplace" "Typical" ...
## $ Garage_Type : chr "Attchd" "Attchd" "Attchd" "Attchd" ...
## $ Garage_Yr_Blt : int 1960 1961 1958 1968 1997 1998 2001 1992 1995 1999 ...
## $ Garage_Finish : chr "Fin" "Unf" "Unf" "Fin" ...
## $ Garage_Cars : int 2 1 1 2 2 2 2 2 2 2 ...
## $ Garage_Area : int 528 730 312 522 482 470 582 506 608 442 ...
## $ Garage_Qual : chr "Typical" "Typical" "Typical" "Typical" ...
## $ Garage_Cond : chr "Typical" "Typical" "Typical" "Typical" ...
## $ Paved_Drive : chr "Partial_Pavement" "Paved" "Paved" "Paved" ...
## $ Wood_Deck_SF : int 210 140 393 0 212 360 0 0 237 140 ...
## $ Open_Porch_SF : int 62 0 36 0 34 36 0 82 152 60 ...
## $ Enclosed_Porch : int 0 0 0 0 0 0 170 0 0 0 ...
## $ Three_season_porch: int 0 0 0 0 0 0 0 0 0 0 ...
## $ Screen_Porch : int 0 120 0 0 0 0 0 144 0 0 ...
## $ Pool_Area : int 0 0 0 0 0 0 0 0 0 0 ...
## $ Pool_QC : chr "No_Pool" "No_Pool" "No_Pool" "No_Pool" ...
## $ Fence : chr "No_Fence" "Minimum_Privacy" "No_Fence" "No_Fence" ...
## $ Misc_Feature : chr "None" "None" "Gar2" "None" ...
## $ Misc_Val : int 0 0 12500 0 0 0 0 0 0 0 ...
## $ Mo_Sold : int 5 6 6 4 3 6 4 1 3 6 ...
## $ Year_Sold : int 2010 2010 2010 2010 2010 2010 2010 2010 2010 2010 ...
## $ Sale_Type : chr "WD " "WD " "WD " "WD " ...
## $ Sale_Condition : chr "Normal" "Normal" "Normal" "Normal" ...
## $ Longitude : num -93.6 -93.6 -93.6 -93.6 -93.6 ...
## $ Latitude : num 42.1 42.1 42.1 42.1 42.1 ...
## $ Sale_Price : int 215000 105000 172000 244000 189900 195500 213500 191500 236500 189000 ...数据集包含2930个observations与83个变量，变量包含数值变量与非数值变量，我们需要在数据预处理时对非数值变量进行处理。
# Error_vec_lasso=rep(0,10)
# Error_vec_xgboost=rep(0,10)
#for (k in 1:10){
k=2
train <- data[-testIDs[,k], ]
test <- data[testIDs[,k], ]
test.y <- test[, c(1, 83)]
test <- test[, -83]
write.csv(train,"train.csv",row.names=FALSE)
write.csv(test, "test.csv",row.names=FALSE)
write.csv(test.y,"test_y.csv",row.names=FALSE)
train <- read.csv("train.csv")
test <- read.csv("test.csv")首先是划分数据集，testIDs中包含打乱的测试数据的序号，用于10-fold检验，在这里将k设置为2，也就是抽取第二组测试数据序号，一共有十组。
train$Garage_Yr_Blt[is.na(train$Garage_Yr_Blt)] = 0
train<- subset(train,select=-c(Street, Utilities,Condition_2, Roof_Matl, Heating, Pool_QC, Misc_Feature, Low_Qual_Fin_SF, Pool_Area, Longitude,Latitude))
test$Garage_Yr_Blt[is.na(test$Garage_Yr_Blt)] = 0
test <- subset(test, select = -c(Street, Utilities,Condition_2, Roof_Matl, Heating, Pool_QC, Misc_Feature, Low_Qual_Fin_SF, Pool_Area, Longitude,Latitude))
test<- test[,-1]首先是查找缺失值，在这个数据集中缺失值只有Garage_Yr_Blt这一列有，由于Garage_Yr_Blt为分类数据，直接将缺失值设置为0，这里的原因是我不想将任何额外信息带入数据集，使用插补法的话，难免会引入原本不存在的信息，但如果将缺失值看作新的一类的话，直观上讲，并没有引入任何信息。
然后我是删除了一些变量，这是因为在直观上，我认为这些变量是与预测变量没有任何关系的，比如说街道，维度，经度以及一些可翻修的房屋变量。
train<- train[,-1]
train.y=train[,dim(train)[2]]
train.x=train[,-dim(train)[2]]
n_train=dim(train.x)[1]
train.x=rbind(train.x,test) ## rbind
categorical.vars <- colnames(train.x)[
which(sapply(train.x,
function(x) mode(x)=="character"))]
train.matrix <- train.x[, !colnames(train.x) %in% categorical.vars,
drop=FALSE]
n.train <- nrow(train.matrix)
for(var in categorical.vars){
mylevels <- sort(unique(train.x[, var]))
m <- length(mylevels)
m <- ifelse(m>2, m, 1)
tmp.train <- matrix(0, n.train, m)
col.names <- NULL
for(j in 1:m){
tmp.train[train.x[, var]==mylevels[j], j] <- 1
col.names <- c(col.names, paste(var, '_', mylevels[j], sep=''))
}
colnames(tmp.train) <- col.names
train.matrix <- cbind(train.matrix, tmp.train)
}
test=train.matrix[c((n_train+1):dim(train.matrix)[1]),]
train.x=train.matrix[c(1:n_train),]然后就是将分类变量变为虚拟变量，也就是变为布尔类型的数据，这无疑会带来维度上的急剧增长，所以对于中大型数据集的虚拟变量的转换，数据降维是必要的。
winsor.vars <- c("Lot_Frontage", "Lot_Area", "Mas_Vnr_Area", "BsmtFin_SF_2", "Bsmt_Unf_SF", "Total_Bsmt_SF", "Second_Flr_SF", 'First_Flr_SF', "Gr_Liv_Area", "Garage_Area", "Wood_Deck_SF", "Open_Porch_SF", "Screen_Porch", "Misc_Val")
quan.value <- 0.95
for(var in winsor.vars){
tmp <- train.x[, var]
myquan <- quantile(tmp, probs = quan.value, na.rm = TRUE)
tmp[tmp > myquan] <- myquan
train.x[, var] <- tmp
}
for(var in winsor.vars){
tmp <- test[, var]
myquan <- quantile(tmp, probs = quan.value, na.rm = TRUE)
tmp[tmp > myquan] <- myquan
test[, var] <- tmp
}下一步就是winsorize数据处理，这是处理离群值点的一种，也是降低过拟合的一种方法。Winsorization是应该将超出变量特定百分位范围的数值替换为其特定百分位数值的方法。
#### scale the numeric variable
train.x$Year_Built=scale(train.x$Year_Built, center = TRUE, scale = FALSE)
train.x$Year_Remod_Add=scale(train.x$Year_Remod_Add, center = TRUE, scale = FALSE)
train.x$Garage_Yr_Blt=scale(train.x$Garage_Yr_Blt, center = TRUE, scale = FALSE)
train.x$Year_Sold=scale(train.x$Year_Sold, center = TRUE, scale = FALSE)
train.x$Year_Built=scale(train.x$Year_Built, center = TRUE, scale = FALSE)
test$Year_Built=scale(test$Year_Built, center = TRUE, scale = FALSE)
test$Year_Remod_Add=scale(test$Year_Remod_Add, center = TRUE, scale = FALSE)
test$Garage_Yr_Blt=scale(test$Garage_Yr_Blt, center = TRUE, scale = FALSE)
test$Year_Sold=scale(test$Year_Sold, center = TRUE, scale = FALSE)
test$Year_Built=scale(test$Year_Built, center = TRUE, scale = FALSE)然后是对某些变量进行归一化，我在这边是取年份变量进行归一化。
par(mfrow=c(2,1))
hist(train.y)
hist(log(train.y))
train.y=log(train.y)上面是我认为做重要的，就是对预测变量取对数。我尝试过很多次，发现这个数据集中，对数化对lasso模型的帮助很大，能大幅提高预测精度，但对于XGBoost的用处几乎没有任何作用。
这里让我们查看一下数据的分布，如上图所示。我查了一些解释，但都没有发现为什么对数化在对于线性回归帮助很大的同时，对于XGBoost等等的树模型的帮助却很小，这个问题我会在后面进行讨论。
下面就是预测模型，用了两个预测模型，一个是通过Lasso进行变量筛选，再进行Ridge回归，一个是XGBoost模型
set.seed(2021)
xgb.model <- xgboost(data = as.matrix(train.x),
label = train.y, max_depth = 6,
eta = 0.05, nrounds = 5000,
subsample = 0.5,
verbose = FALSE)
cv.out <- cv.glmnet(as.matrix(train.x[, xgb.model$feature_names]), train.y, alpha = 1)
sel.vars <- predict(cv.out, type="nonzero",
s = cv.out$lambda.1se)[,1]
cv.out <- cv.glmnet(as.matrix(train.x[, sel.vars]),
train.y, alpha = 0)
test.y <- read.csv("test_y.csv")
test.y=test.y[,2]
Ytest.pred <-predict(cv.out, s = cv.out$lambda.min, newx = as.matrix(test[, sel.vars]))
Ytest.pred=exp(Ytest.pred)
Error_vec_lasso=sqrt(mean((log(test.y) - log(Ytest.pred))^2))
xgb.pred = predict(xgb.model, newdata = as.matrix(test[,xgb.model$feature_names]))
xgb.pred=exp(xgb.pred)
Error_vec_xgboost=sqrt(mean((log(test.y) - log(xgb.pred))^2))
#}
Error_vec_lasso## [1] 0.120437Error_vec_xgboost## [1] 0.1194583得到了两个结果，都在0.125以下，比较下来结果比较满意。写一个for循环就可以得到10组两个模型的结果，这里就不展示了