-
Notifications
You must be signed in to change notification settings - Fork 3
/
Copy pathLinearRegression.R
143 lines (103 loc) · 5.19 KB
/
LinearRegression.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
#Linear regression Step by Step
#-------------------------------------------------------------------#
#library
#-------------------------------------------------------------------#
library(car)
library(lmtest)
library(ggplot2)
#-------------------------------------------------------------------#
#input parameters
#-------------------------------------------------------------------#
#path pointing to the input .csv file
c_path_in <- "D:/Dataset/sample.csv"
#output folder path
c_path_out <- "D:/temp/"
#continuous independent variables
c_var_in_independent <- c('var_1','var_2','var_3','var_4')
#continuous dependent variable
c_var_in_dependent <- c('var_5')
#-------------------------------------------------------------------#
#function
#-------------------------------------------------------------------#
mape <- function(actual,predicted) mean(abs((actual - predicted)/actual))*100
rmse <- function (actual, predicted) sqrt(mean(actual-predicted)^2)
#-------------------------------------------------------------------#
#load input dataset
#-------------------------------------------------------------------#
data <- read.csv(c_path_in)
#-------------------------------------------------------------------#
#subset
#-------------------------------------------------------------------#
data <- data[,c(c_var_in_independent,c_var_in_dependent)]
data <- na.omit(data)
#-------------------------------------------------------------------#
#model building
#-------------------------------------------------------------------#
formula <- as.formula(paste0(c_var_in_dependent, '~' , paste0(c_var_in_independent,collapse = '+')))
lmObj <- lm(formula = formula,
data = data)
#-------------------------------------------------------------------#
#model information
#-------------------------------------------------------------------#
parameterEstimate <- coef(summary(lmObj))
varianceInflationFactor <- car::vif(lmObj)
varianceInflationFactor <- c(NA,varianceInflationFactor)
actual <- lmObj$model[,c_var_in_dependent]
predicted <- lmObj$fitted.values
residual <- round(residuals(lmObj),4)
mape <- mape(actual = actual,
predicted = predicted)
RSquared <- summary(lmObj)$r.squared
adjustedRSquared <- summary(lmObj)$adj.r.squared
FStatistic <- summary(lmObj)$fstatistic["value"]
rootMeanSquaredError <- rmse(actual, predicted)
#-------------------------------------------------------------------#
#Goldfeld Quandt test for Homoscedasticity
#-------------------------------------------------------------------#
goldfledQuandtTest <- lmtest::gqtest(lmObj)[1]
#-------------------------------------------------------------------#
#Durbin watson test for autocorrelation (1.5 < stat < 2.5)
#-------------------------------------------------------------------#
durbinWatsonTest <- lmtest::dwtest(lmObj)[1]
#-------------------------------------------------------------------#
#ouput
#-------------------------------------------------------------------#
parameterEstimate <- cbind.data.frame(parameterEstimate,varianceInflationFactor)
modelStatistic <- cbind.data.frame(Statistic = c("R-Squared","Adj. R-Squared","Root Mean Squared Error","Mean Absolute Percentage Error","F-Statistic"),
Value =c(RSquared,adjustedRSquared,rootMeanSquaredError,mape,FStatistic))
write.csv(parameterEstimate,paste0(c_path_out,"/Estimates.csv"))
write.csv(modelStatistic,paste0(c_path_out,"/ModelSatistic.csv"))
#-------------------------------------------------------------------#
#Actual vs Predicted
#-------------------------------------------------------------------#
png(filename=paste0(c_path_out,"/ActualPredicted.png"),
width = 12,
height = 6,
units = "in",
res =100
)
rs <- qplot(x = actual,
y = predicted,
xlab = "Actual",
ylab = "Predicted") + theme(aspect.ratio = 1/2.5)
print({
rs2 <- rs + theme(axis.text.x=element_text(size=12, angle=40, vjust=.9, hjust=1.01))
})
dev.off()
#-------------------------------------------------------------------#
#Residual vs Predicted plot verifying the assumptions of Linear Model
#-------------------------------------------------------------------#
png(filename = paste0(c_path_out,"/ResidualPredicted.png"),
width = 12,
height = 6,
units = "in",
res = 100
)
rs <- qplot(x = predicted,
y = residual,
xlab = "Predicted",
ylab = "Residual") + theme(aspect.ratio = 1/2.5)
print({
rs2 <- rs + theme(axis.text.x=element_text(size=12, angle=40, vjust=.9, hjust=1.01))
})
dev.off()