R Exam Notebook: MN50645

University of Bath, School Of Management, MN50645 Data Mining 2019/20

Student id: 199042788

Name: Tianyun Cai

---

Business Analytics Data Case Study: Property Market

Q1 (a)

---

Downloading and preparing data

#Load packages
library(xray)
library(corrplot)
library(mice)
library(sqldf)
library(lattice)
library(MASS)
library(nnet)
library(ggplot2)
library(ggthemes)
library(neuralnet)
library(NeuralNetTools)

#Set working directory
setwd("C:/Users/skych/OneDrive - University of Bath/Desktop/COURSE MATERIAL/MN50645 Data mining/assessment")

#Read file from working directory
property <- read.csv("MN50645_exam_dataset[23].csv", header=T)

---

Firstly, we apply str() to see the structure of this dataset. As shown below, there are 10000 observations of 13 variables.

str(property)

'data.frame':   10000 obs. of  13 variables:
 $ Suburb       : Factor w/ 327 levels "Abbotsford","Aberfeldie",..: 275 132 144 53 256 156 192 12 107 155 ...
 $ Rooms        : int  3 3 4 4 4 3 3 4 4 6 ...
 $ Type         : Factor w/ 3 levels "h","t","u": 2 1 1 1 1 1 1 1 1 1 ...
 $ Price        : int  520000 NA NA 1225000 920000 NA 695000 2000000 854000 2825000 ...
 $ Method       : Factor w/ 9 levels "PI","PN","S",..: 1 3 3 3 6 5 3 1 3 3 ...
 $ Distance     : num  14 3.5 16.7 11.8 23.5 7.5 12.7 6.3 22.2 5.3 ...
 $ Bedroom2     : int  3 3 4 4 4 3 NA 4 4 6 ...
 $ Bathroom     : int  2 1 2 2 2 2 NA 3 2 5 ...
 $ Car          : int  2 3 2 2 2 2 NA 4 2 2 ...
 $ Landsize     : int  240 382 689 723 727 229 NA 544 676 NA ...
 $ BuildingArea : num  134 103 NA 159 218 NA NA NA 181 NA ...
 $ YearBuilt    : int  2007 1900 NA 1970 1990 1900 NA NA 1975 1920 ...
 $ Propertycount: int  14042 6244 15321 4480 11667 6482 4168 4836 3940 11308 ...

---

Then, the function of anomalies() can help us detection of anomalies and unusual values. The variables of “BuildingArea” and “YearBuilt” have missing values exceeding 50%; The variable of “Car” and “Bathroom” have 33.59% and 24.51% missing values and over 4% zore values; And the variable of “Bedroom2” and “Price” have over 20% missing values. In general, there is a great number of null values and some zero values in this data set that need further processing.

anomalies(property)

$variables

$problem_variables

NA

---

Q1 (b)

---

Before fixing the problem of our dataset, we should determine the importance of the problematic variables. As we should explorate the factors affecting propert prices and also build a ML forecasting model to predict price in the Q2 and Q3, we will caculate the corelartion bewteen the problematic variables and the price first.

property_cor1 <- cor(na.omit(sqldf("select Price, BuildingArea from property")))
property_cor2 <- cor(na.omit(sqldf("select Price, Landsize from property")))
property_cor3 <- cor(na.omit(sqldf("select Price, YearBuilt from property")))
property_cor4 <- cor(na.omit(sqldf("select Price, Car from property")))
property_cor5 <- cor(na.omit(sqldf("select Price, Bathroom from property")))
property_cor6 <- cor(na.omit(sqldf("select Price, Bedroom2 from property")))
par(mfcol=c(2,3))
corrplot.mixed(property_cor1, lower = "number", upper = "pie", tl.col = "black",lower.col = "black", number.cex = 1)
corrplot.mixed(property_cor2, lower = "number", upper = "pie", tl.col = "black",lower.col = "black", number.cex = 1)

corrplot.mixed(property_cor3, lower = "number", upper = "pie", tl.col = "black",lower.col = "black", number.cex = 1)
corrplot.mixed(property_cor4, lower = "number", upper = "pie", tl.col = "black",lower.col = "black", number.cex = 1)

corrplot.mixed(property_cor5, lower = "number", upper = "pie", tl.col = "black",lower.col = "black", number.cex = 1)
corrplot.mixed(property_cor6, lower = "number", upper = "pie", tl.col = "black",lower.col = "black", number.cex = 1)

---

As shown in above figure, “Landsize” and “Car” has few relevant to price while it has over 20% missing value, so we will give them up. For the rest problematic variables, we consider using appropriate methods to fill in missing values.

property <- property[,c(-1, -9,-10)]

---

For the missing value in “Price”, we will delete the missing value because “Price” is our dependent variable. And for “Propertycount”, we will also delete its missing value since the few number.

mvr1 = which(is.na(property$Price))
property <- property[-mvr1,]
mvr2 = which(is.na(property$Propertycount))
property <- property[-mvr2,]

---

By visualizing “YearBuilt”, we can see that most properties were built bewtween 1950 and 2010. Its mode is in 2000-2010. However, its median and mean value is about 1970. Thus, median may be more representative.

ggplot(property, aes(x = YearBuilt)) + geom_histogram(bins = 20, fill = "lightblue", colour = "black")

summary(property$YearBuilt)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
   1830    1950    1970    1968    2000    2018    4272 

---

Therefore, we will fill with median.

property[is.na(property$YearBuilt),"YearBuilt"] <- median(property$YearBuilt,na.rm = T)

---

Factors such as the number of bedrooms, bathrooms, house size, room type, and price could be all related to each other. Thus, for “BuildingArea”, “Bathroom” and “Bedroom2”, we will delete the column with missing value of these three variables first.

property <-sqldf("select * from property 
              where BuildingArea is not null
              or Bathroom is not null
              or Bedroom2 is not null",row.names=TRUE)

---

Then we check the anomalies and unusual values again.

anomalies(property)

$variables

$problem_variables

NA

---

There is only 1 missing value in “Bathroom” which we can just delete and the missing value in “BuildingArea” is still occupy 48.69%. Hence, we should filling function to deal with the problem. Besides, we can aslo find there are 16 zero value in “BuildingArea”, which is unreasonable, we will also delete the line with the variable less than 10.

mvr3 = which(is.na(property$Bathroom))
property <- property[-mvr3,]
zvr1 = which(property$BuildingArea < 10)
property <- property[-zvr1,]

---

Mice() function can help use reasonable data values to estimate missing values. We try use all avriables to build model first.

imp<-mice(property, method="rf", seed=1234)  
fit<-with(imp,lm(Price~.,data=property)) 
pool=pool(fit) 
options(digits=3)  
summary(pool) 

---

As shown above, the p-values of “Method”, “Propertycount” and “Bedroom2” are not significant. So we will delete them and try again. The result shows that this model could do a good job.

imp<-mice(property[,c(1,2,3,5,7,8,9)], method="rf", seed=1234)  

 iter imp variable
  1   1  BuildingArea
  1   2  BuildingArea
  1   3  BuildingArea
  1   4  BuildingArea
  1   5  BuildingArea
  2   1  BuildingArea
  2   2  BuildingArea
  2   3  BuildingArea
  2   4  BuildingArea
  2   5  BuildingArea
  3   1  BuildingArea
  3   2  BuildingArea
  3   3  BuildingArea
  3   4  BuildingArea
  3   5  BuildingArea
  4   1  BuildingArea
  4   2  BuildingArea
  4   3  BuildingArea
  4   4  BuildingArea
  4   5  BuildingArea
  5   1  BuildingArea
  5   2  BuildingArea
  5   3  BuildingArea
  5   4  BuildingArea
  5   5  BuildingArea

fit<-with(imp,lm(Price~.,data=property[,c(1,2,3,5,7,8,9)])) 
pool=pool(fit) 
options(digits=3)  
summary(pool)

          term estimate std.error statistic   df  p.value
1  (Intercept) 10108105    543209     18.61 3044 0.00e+00
2        Rooms    84002     13374      6.28 3044 3.85e-10
3        Typet   -97499     31985     -3.05 3044 2.32e-03
4        Typeu  -203481     28782     -7.07 3044 1.92e-12
5     Distance   -29197      1347    -21.68 3044 0.00e+00
6     Bathroom   235531     15961     14.76 3044 0.00e+00
7 BuildingArea     2512       130     19.37 3044 0.00e+00
8    YearBuilt    -4925       283    -17.43 3044 0.00e+00

property_mice=complete(imp,action=1) 

---

Finally, we update the property data and check the anomalies and unusual values again. There is no missing value remaining and only two variables have zero values.

property[,8] <- property_mice[,6]

anomalies(property)

$variables

$problem_variables

NA

---

Q2 (a)

---

Basic descriptive statistics with summary() function. From the result, we can get a brief understand of the dataset. “Type” and “Method” are string format, which should do some transform later.

summary(property)

     Rooms       Type         Price              Method        Distance       Bedroom2    
 Min.   : 1.00   h:4511   Min.   :  145000   S      :3917   Min.   : 0.0   Min.   : 0.00  
 1st Qu.: 2.00   t: 475   1st Qu.:  654500   SP     : 744   1st Qu.: 6.4   1st Qu.: 2.00  
 Median : 3.00   u: 998   Median :  903000   PI     : 718   Median :10.5   Median : 3.00  
 Mean   : 3.07            Mean   : 1083052   VB     : 558   Mean   :11.4   Mean   : 3.05  
 3rd Qu.: 4.00            3rd Qu.: 1330875   SA     :  47   3rd Qu.:14.2   3rd Qu.: 4.00  
 Max.   :12.00            Max.   :11200000   PN     :   0   Max.   :47.3   Max.   :12.00  
                                             (Other):   0                                 
    Bathroom    BuildingArea    YearBuilt    Propertycount  
 Min.   :0.0   Min.   :  11   Min.   :1830   Min.   :   83  
 1st Qu.:1.0   1st Qu.: 101   1st Qu.:1966   1st Qu.: 4380  
 Median :1.5   Median : 130   Median :1970   Median : 6567  
 Mean   :1.6   Mean   : 145   Mean   :1969   Mean   : 7552  
 3rd Qu.:2.0   3rd Qu.: 173   3rd Qu.:1977   3rd Qu.:10412  
 Max.   :8.0   Max.   :1561   Max.   :2018   Max.   :21650  
                                                            

---

To make above result more easy to unstanding, we can visualize them. The figure below shows the density distribution of BuildingArea variable. Most properties’s building size are bewteen 100 and 140. Only few of them larger than 400.

property_temp1 <- filter(property, BuildingArea  <= 600)
ggplot(property_temp1, aes(x = BuildingArea)) +
  geom_histogram(aes(y = ..density..), fill = "lightblue", binwidth = 20, color="gray") + 
  geom_density() +
  theme_pander()+
  ggtitle("BuildingArea Distribution")+
  theme(plot.title = element_text(hjust = 0.5))

---

For the methods of the property sold, “property sold” occupies the majority which is almost 4000. “property sold prior” and “property passed in” have about 700, while “vendor bid” only have less than 600. For types of the preporty, h (house, cottage, villa, semi, terrace) occupies the majority. u (unit, duplex) hasabout 1000 and townhouse just have about 500.

ggplot(property) +
  geom_bar(aes(x = reorder(Method, table(Method)[Method])), fill = "lightblue") +
  coord_flip() +
  theme_bw() +
  labs(x = NULL)

ggplot(property) +
  geom_bar(aes(x = reorder(Type, table(Type)[Type])), fill = "lightblue", width = 0.7) +
  coord_flip() +
  theme_bw() +
  labs(x = NULL)

---

Then we plot the density distribution of price of different methods and types. We can see that different “Method” and “Type” have the almost similiar mode and median, which means “Method” and “Type” may have a few influence of price.

ggplot(property, aes(x = Price)) + 
  geom_density(fill="yellow") +
  facet_wrap( ~ Method) +
  geom_vline(aes(xintercept = median(Price)), linetype = 2) +
  theme_bw()

ggplot(property, aes(x = Price)) + 
  geom_density(fill="yellow") +
  facet_wrap( ~ Type) +
  geom_vline(aes(xintercept = median(Price)), linetype = 2) +
  theme_bw()

---

Q2 (b)

---

Question: 1. How distance influence the sales of property? 2. How “Method” influence the sales of property? 3. How “Type” influence the sales of property?

---

To solve the question 1, we plot the point plot to show the relationship between “Distance” and “Price” as well as a fitted line. It shows that distance will influence the price of the property. So the closer the property is to the city center, the higher the demand is (According to the Economics, demand determines price).

ggplot(property, aes(x = Distance, y = Price)) + 
  geom_point(size = 0.5, color="darkgray") +
  geom_smooth(method = "lm") +
  theme_bw()

---

As mentioned above, “BuildingArea” has the highest correlation with price. Then we could explore how building size influence the price of different “Method” and “Type”. From the fitted curve in two figures below, we can find that sold after auction and vendor bid properties’ building size have a higher unit price (price per size). For different types, townhouse may have a higher unit price. Thus, these features may improve the sales of property when the total area of the house sold is certain.

ggplot(property, aes(x = BuildingArea, y = Price)) + 
  geom_point(size = 0.5, color="darkgray") +
  geom_smooth(method = "lm") +
  facet_wrap( ~ Method) +
  theme_bw()

ggplot(property, aes(x = BuildingArea, y = Price)) + 
  geom_point(size = 0.5, color="darkgray") +
  geom_smooth(method = "lm") +
  facet_wrap( ~ Type) +
  theme_bw()

---

Q3 (a)

---

In this part, we need to build a machine learning model to predict the factors that affect the sales price. In machine learning models that solve numerical problems, Artificial Neural Network algorithms may have better performance. Its classification accuracy is high and its learning ability is extremely strong. It is robust and fault-tolerant to noisy data. In addition, it has the ability to associate and can approximate any nonlinear relationship. ***

Q3 (b)

---

As metioned above, we remoted the two variables “Lansdsize” and “Car” as they are not highly correlated with “Prices” and have more vacancies. Then we apply ggpairs() function to check the correlatin among these variable. The correlation between “Propertycount” and “Price” is only about 5%. The number of properties that exist in the suburb actually can not influent the property price in real world. To make our model more concise and efficient, we will delete this variable.

ggpairs(property, title = "Scatterplot Matrix of the Features of the Property Data Set")

property <- property[,-10]

---

Q3 (c)

---

Before applying the ANN function, we should transform “Method” and “Type” into numeric format.

property_ANN <- property

t = array(0, 5984)
h = array(0, 5984)
u = array(0, 5984)
S = array(0, 5984)
SP = array(0, 5984)
PI = array(0, 5984)
VB = array(0, 5984)
SA = array(0, 5984)

for (i in 1: 5984) {
  
  if(property_ANN$Type[i] == "t") {t[i]=1  
  h[i]=0
  u[i]=0}  
  
  if(property_ANN$Type[i] == "h") {t[i]=0
  h[i]=1
  u[i]=0}  
  
  if(property_ANN$Type[i] == "u") {t[i]=0
  h[i]=0
  u[i]=1}
  
  if(property_ANN$Method[i] == "S") {S[i]=1  
  SP[i]=0
  PI[i]=0
  VB[i]=0
  SA[i]=0}
  
  if(property_ANN$Method[i] == "SP") {S[i]=0  
  SP[i]=1
  PI[i]=0
  VB[i]=0
  SA[i]=0}
  
  if(property_ANN$Method[i] == "PI") {S[i]=0  
  SP[i]=0
  PI[i]=1
  VB[i]=0
  SA[i]=0} 
  
  if(property_ANN$Method[i] == "VB") {S[i]=0  
  SP[i]=0
  PI[i]=0
  VB[i]=1
  SA[i]=0} 
  
  if(property_ANN$Method[i] == "SA") {S[i]=0  
  SP[i]=0
  PI[i]=0
  VB[i]=0
  SA[i]=1} 
}
property_ANN$Type_t=t
property_ANN$Type_h=h
property_ANN$Type_u=u
property_ANN$Method_S=S
property_ANN$Method_SP=SP
property_ANN$Method_PI=PI
property_ANN$Method_VB=VB
property_ANN$Method_SA=SA
property_ANN <- property_ANN[,c(-2,-4)]

---

Data observation should be scaled onto the [0,1] interval.

scale01 <- function(x){
  (x - min(x)) / (max(x) - min(x))
}
property_ANN <- property_ANN %>%mutate_all(scale01)

---

Then spliting into test(0.2) and train sets (0.8).

set.seed(12341)
Property_Train <- sample_frac(tbl = property_ANN, replace = FALSE, size = 0.80)
Property_Test <- anti_join(property_ANN, Property_Train)

---

Firstly, we construct a model with 1-hidden layer. By using plot() function, we can get the weights learned by the ANN_model1 neural network.The plot also shows the sum of Squared Errors (SSE) of the training data set and displays the number of training iterations before convergence. In this model, SSE of training data set is 0.7018.

set.seed(12342)
ANN_model1 <- neuralnet(Price ~ ., data = Property_Train)
plot(ANN_model1)

---

Then we visualise the weights of variables. We can see that the “BuildingArea” has the highest importance of the price pridicting, so it may influence price most. “Distance” and “Bathroom” is the second and third importance factor. “Rooms” and “Type_u” have over 0.05 importance.

garson(ANN_model1) + theme(axis.text.x = element_text(angle = 45, hjust = 1))

plotnet(ANN_model1)

---

Q3 (d)

---

To evaluate the performance, we will use SSE, MSE and R-squared.

#create function to caculate performance metrics
ssefun <- function(pred, obs) sum((as.numeric(unlist(pred)) - obs)^2)/2
msefun <- function(pred, obs) mean((as.numeric(unlist(pred)) - obs)^2/2)
sstfun <- function(obs) sum((obs - mean(obs))^2)/2

---

By caculating the three metrics of the training model. According to the result, we can find that the SSE is 4.1184, Mse is 0.0009 and R-squared is 0.4948, which suggests results can be explained by the train model. Since R-squared is not so close to 1, it should get improved.

sse11 <- ssefun(ANN_model1$net.result, Property_Train[, 2])
mse11 <- msefun(ANN_model1$net.result, Property_Train[, 2])
rs11 <- 1-sse11/sstfun(Property_Train[, 2])
paste("Train SSE: ", round(sse11, 4), "; Train MSE: ", round(mse11, 4), "; Train R-squared: ", round(rs11, 4))

[1] "Train SSE:  4.1184 ; Train MSE:  9e-04 ; Train R-squared:  0.4948"

---

Then we run test observations and caculating the performance metrics. The SSE of test model is better than train model, which is only 1.0555. However, MSE did not get a significant improvement, so actually the test model has almost the same performance with the training one, which can also be found in R-squared.

#Run test observations
Test_ANN1_Output <- compute(ANN_model1, Property_Test[, -2])$net.result

sse12 <- ssefun(Test_ANN1_Output, Property_Test[, 2])
mse12 <- msefun(Test_ANN1_Output, Property_Test[, 2])
rs12 <- 1-sse12/sstfun(Property_Test[, 2])
paste("Test SSE: ", round(sse12, 4), "; Test MSE: ", round(mse12, 4), "; Test R-squared: ", round(rs12, 4))

[1] "Test SSE:  1.0555 ; Test MSE:  9e-04 ; Test R-squared:  0.5231"

---

Q3 (e)

---

To improve our model, we could use “rep” parameter to let the ANN algrithnm run many time and find the best one.

set.seed(12343)
ANN_model2 <- neuralnet(Price ~ ., data = Property_Train, rep = 5)

Warning in for (w in (length.weights - 1):1) { :
  closing unused connection 5 (<-DESKTOP-KMTTC50:11822)
Warning in for (w in (length.weights - 1):1) { :
  closing unused connection 4 (<-DESKTOP-KMTTC50:11822)
Warning in for (w in (length.weights - 1):1) { :
  closing unused connection 3 (<-DESKTOP-KMTTC50:11822)
Algorithm did not converge in 1 of 5 repetition(s) within the stepmax.

plot(ANN_model2, rep = "best")

Min_num2 = which.min(ANN_model2$result.matrix[1,])

---

Then we caculate the performance of this train model and test model.

sse21 <- ssefun(ANN_model2$net.result[[Min_num2]], Property_Train[, 2])
mse21 <- msefun(ANN_model2$net.result[[Min_num2]], Property_Train[, 2])
rs21 <- 1-sse21/sstfun(Property_Train[, 2])
paste("Train SSE: ", round(sse21, 4), "; Train MSE: ", round(mse21, 4), "; Train R-squared: ", round(rs21, 4))

[1] "Train SSE:  4.1149 ; Train MSE:  9e-04 ; Train R-squared:  0.4952"

Test_ANN2_Output <- compute(ANN_model2, Property_Test[, -2], rep = Min_num2)$net.result

sse22 <- ssefun(Test_ANN2_Output, Property_Test[, 2])
mse22 <- msefun(Test_ANN2_Output, Property_Test[, 2])
rs22 <- 1-sse22/sstfun(Property_Test[, 2])
paste("Test SSE: ", round(sse22, 4), "; Test MSE: ", round(mse22, 4), "; Test R-squared: ", round(rs22, 4))

[1] "Test SSE:  1.0474 ; Test MSE:  9e-04 ; Test R-squared:  0.5268"

---

In addition, we can also increase the hidden neurons number into 10 (15 variables * 2/3). Both train model and test model get improved. 10 hidden could do a better job than the former one.

set.seed(12344)
ANN_model3 <- neuralnet(Price ~ ., data = Property_Train, hidden = 10)

sse31 <- ssefun(ANN_model3$net.result, Property_Train[, 2])
mse31 <- msefun(ANN_model3$net.result, Property_Train[, 2])
rs31 <- 1-sse31/sstfun(Property_Train[, 2])
paste("Train SSE: ", round(sse31, 4), "; Train MSE: ", round(mse31, 4), "; Train R-squared: ", round(rs31, 4))

[1] "Train SSE:  3.6439 ; Train MSE:  8e-04 ; Train R-squared:  0.553"

Test_ANN3_Output <- compute(ANN_model3, Property_Test[, -2])$net.result

sse32 <- ssefun(Test_ANN3_Output, Property_Test[, 2])
mse32 <- msefun(Test_ANN3_Output, Property_Test[, 2])
rs32 <- 1-sse32/sstfun(Property_Test[, 2])
paste("Test SSE: ", round(sse32, 4), "; Test MSE: ", round(mse32, 4), "; Test R-squared: ", round(rs32, 4))

[1] "Test SSE:  0.9688 ; Test MSE:  8e-04 ; Test R-squared:  0.5622"

---

We can also add layer into 2 (2-Hidden Layers, Layer-1 10-neurons, Layer-2, 3-neuron). The train model is better than the former one, but test model does not improve.

set.seed(12345)
ANN_model4 <- neuralnet(Price ~ ., data = Property_Train, hidden = c(10, 3), stepmax = 10000000)

sse41 <- ssefun(ANN_model4$net.result, Property_Train[, 2])
mse41 <- msefun(ANN_model4$net.result, Property_Train[, 2])
rs41 <- 1-sse41/sstfun(Property_Train[, 2])
paste("Train SSE: ", round(sse41, 4), "; Train MSE: ", round(mse41, 4), "; Train R-squared: ", round(rs41, 4))

[1] "Train SSE:  3.2109 ; Train MSE:  7e-04 ; Train R-squared:  0.6061"

Test_ANN4_Output <- compute(ANN_model4, Property_Test[, -2])$net.result

sse42 <- ssefun(Test_ANN4_Output, Property_Test[, 2])
mse42 <- msefun(Test_ANN4_Output, Property_Test[, 2])
rs42 <- 1-sse42/sstfun(Property_Test[, 2])
paste("Test SSE: ", round(sse42, 4), "; Test MSE: ", round(mse42, 4), "; Test R-squared: ", round(rs42, 4))

[1] "Test SSE:  1.0043 ; Test MSE:  8e-04 ; Test R-squared:  0.5462"

---

Changing activation function could be a way to improve our model. This model has higher R-squared and lower SSE/MSE in test model, but the improvement is not significant.

set.seed(12346)
ANN_model5 <- neuralnet(Price ~ ., data = Property_Train, hidden = 10, act.fct = "logistic")

sse51 <- ssefun(ANN_model5$net.result, Property_Train[, 2])
mse51 <- msefun(ANN_model5$net.result, Property_Train[, 2])
rs51 <- 1-sse51/sstfun(Property_Train[, 2])
paste("Train SSE: ", round(sse51, 4), "; Train MSE: ", round(mse51, 4), "; Train R-squared: ", round(rs51, 4))

[1] "Train SSE:  3.5889 ; Train MSE:  7e-04 ; Train R-squared:  0.5597"

Test_ANN5_Output <- compute(ANN_model5, Property_Test[, -2])$net.result

sse52 <- ssefun(Test_ANN5_Output, Property_Test[, 2])
mse52 <- msefun(Test_ANN5_Output, Property_Test[, 2])
rs52 <- 1-sse52/sstfun(Property_Test[, 2])
paste("Test SSE: ", round(sse52, 4), "; Test MSE: ", round(mse52, 4), "; Test R-squared: ", round(rs52, 4))

[1] "Test SSE:  0.9418 ; Test MSE:  8e-04 ; Test R-squared:  0.5745"

---

Therefore, we will try more models to find the best one.

set.seed(12347)
#hidden = 8
ANN_model6 <- neuralnet(Price ~ ., data = Property_Train, hidden = 8)
sse61 <- ssefun(ANN_model6$net.result, Property_Train[, 2])
mse61 <- msefun(ANN_model6$net.result, Property_Train[, 2])
rs61 <- 1-sse61/sstfun(Property_Train[, 2])
Test_ANN6_Output <- compute(ANN_model6, Property_Test[, -2])$net.result
sse62 <- ssefun(Test_ANN6_Output, Property_Test[, 2])
mse62 <- msefun(Test_ANN6_Output, Property_Test[, 2])
rs62 <- 1-sse62/sstfun(Property_Test[, 2])

#hidden = c(8, 3)
ANN_model7 <- neuralnet(Price ~ ., data = Property_Train, hidden = c(8, 3))
sse71 <- ssefun(ANN_model7$net.result, Property_Train[, 2])
mse71 <- msefun(ANN_model7$net.result, Property_Train[, 2])
rs71 <- 1-sse71/sstfun(Property_Train[, 2])
Test_ANN7_Output <- compute(ANN_model7, Property_Test[, -2])$net.result
sse72 <- ssefun(Test_ANN7_Output, Property_Test[, 2])
mse72 <- msefun(Test_ANN7_Output, Property_Test[, 2])
rs72 <- 1-sse72/sstfun(Property_Test[, 2])

#hidden = c(10, 3), act.fct = "logistic"
ANN_model8 <- neuralnet(Price ~ ., data = Property_Train, hidden = c(10, 3), act.fct = "logistic")
sse81 <- ssefun(ANN_model8$net.result, Property_Train[, 2])
mse81 <- msefun(ANN_model8$net.result, Property_Train[, 2])
rs81 <- 1-sse81/sstfun(Property_Train[, 2])
Test_ANN8_Output <- compute(ANN_model8, Property_Test[, -2])$net.result
sse82 <- ssefun(Test_ANN8_Output, Property_Test[, 2])
mse82 <- msefun(Test_ANN8_Output, Property_Test[, 2])
rs82 <- 1-sse82/sstfun(Property_Test[, 2])

#Rescale for tanh activation function
scale11 <- function(x) {
  (2 * ((x - min(x))/(max(x) - min(x)))) - 1
}
Property_Train1 <- Property_Train %>% mutate_all(scale11)
Property_Test1 <- Property_Test %>% mutate_all(scale11)

#hidden = 10, act.fct = "tanh"
ANN_model9 <- neuralnet(Price ~ ., data = Property_Train1, act.fct = "tanh")
sse91 <- ssefun(ANN_model9$net.result, Property_Train1[, 2])
mse91 <- msefun(ANN_model9$net.result, Property_Train1[, 2])
rs91 <- 1-sse91/sstfun(Property_Train1[, 2])
Test_ANN9_Output <- compute(ANN_model9, Property_Test1[, -2])$net.result
sse92 <- ssefun(Test_ANN9_Output, Property_Test1[, 2])
mse92 <- msefun(Test_ANN9_Output, Property_Test1[, 2])
rs92 <- 1-sse92/sstfun(Property_Test1[, 2])

---

Then we visualize the metrics of these models. According to the results, ANN4 model perform better in train model while ANN8 model do a better job in test model.

tibble(Network = rep(c("ANN1", "ANN2", "ANN3", "ANN4", "ANN5", "ANN6", "ANN7", "ANN8",
                       "ANN9"), each = 2), 
       DataSet = rep(c("Train SSE", "Test SSE"), time = 9), 
       Metrics = c(sse11, sse12, sse21, sse22, sse31, sse32, sse41, sse42, sse51, sse52,
                   sse61, sse62, sse71, sse72, sse81, sse82, sse91, sse92)) %>%
  ggplot(aes(Network, Metrics, fill = DataSet)) + 
  geom_text(aes(label=round(Metrics,3)), position=position_dodge(0.9), vjust = -0.9) + 
  geom_col(position = "dodge") + theme_bw() + ggtitle("Distribution of ANN's SSE")

tibble(Network = rep(c("ANN1", "ANN2", "ANN3", "ANN4", "ANN5", "ANN6", "ANN7", "ANN8",
                       "ANN9"), each = 2), 
       DataSet = rep(c("Train MSE", "Test MSE"),
                     time = 9), 
       Metrics = c(mse11, mse12, mse21, mse22, mse31, mse32, mse41, mse42, mse51, mse52,
                   mse61, mse62, mse71, mse72, mse81, mse82, mse91, mse92)) %>%
  ggplot(aes(Network, Metrics, fill = DataSet)) + 
  geom_text(aes(label=round(Metrics,3)), position=position_dodge(0.9), vjust = -0.9) + 
  geom_col(position = "dodge") + theme_bw() + ggtitle("Distribution of ANN's MSE")

tibble(Network = rep(c("ANN1", "ANN2", "ANN3", "ANN4", "ANN5", "ANN6", "ANN7", "ANN8",
                       "ANN9"), each = 2), 
       DataSet = rep(c("Train R2", "Test R2"),
                     time = 9), 
       Metrics = c(rs11, rs12, rs21, rs22, rs31, rs32, rs41, rs42, rs51, rs52, rs61, rs62,
                   rs71, rs72, rs81, rs82, rs91, rs92)) %>%
  ggplot(aes(Network, Metrics, fill = DataSet)) + 
  geom_text(aes(label=round(Metrics,3)), position=position_dodge(0.9), vjust = -0.9) + 
  geom_col(position = "dodge") + theme_bw() + ggtitle("Distribution of ANN's R-squared")

---

We plot the relationship map of ANN4.

plot(ANN_model4)

---

We plot the relationship map of ANN8.

plot(ANN_model8)

---

Due to the limited number of models we try, this method is not optimal. We can further try deep learning, such as convolutional neural networks. On the basis of the original multi-layer neural network, a feature learning part is added, which imitates the classification of signal processing by the human brain.

---

Q3 (f)

---

According to the results of data exploration (Figure 1 and Figure 2) and meachine learning model(Figure 3), “BuildingArea”, “Bathroom” and “Distance” have highly influence of the property price. So if the agency want to get more sales amount, they should consider get more property with bigger building size, more bathrooms and less distance to city center. And if they want to get more return of invest (ROI), unit price should be a important metrics. Thus, they should use more method of sold after auction and vendor bid. Townhouse will also help them to achieve it.

#Figure 1
ggplot(property, aes(x = BuildingArea, y = Price)) + 
  geom_point(size = 0.5, color="darkgray") +
  geom_smooth(method = "lm") +
  facet_wrap( ~ Method) +
  theme_bw()

ggplot(property, aes(x = BuildingArea, y = Price)) + 
  geom_point(size = 0.5, color="darkgray") +
  geom_smooth(method = "lm") +
  facet_wrap( ~ Type) +
  theme_bw()

#Figure 2
ggplot(property, aes(x = Distance, y = Price)) + 
  geom_point(size = 0.5, color="darkgray") +
  geom_smooth(method = "lm") +
  theme_bw()

#Figure 3
garson(ANN_model1) + theme(axis.text.x = element_text(angle = 45, hjust = 1))