library(MASS)
library(ggplot2)
library(data.table)



#path="C:/Users/user/Dropbox/Michi_Zhenia/Move_to_Atlanta/in_Atlanta/Teaching_at_GaTech/Econometrics_1/lectures_mkm/lectures_computer/computer_session_1/data/"
path_to_folder="C:/User/Kummer/Dropbox/michi_zhenia/Norwich/UEA/Teaching_UEA/ECOnometrics-R002/Lecture_2020_UEA_Short/Computer"
path_to_folder="C:\Users\Kummer\Dropbox\michi_zhenia\Norwich\UEA\Teaching_UEA\ECOnometrics-R002\Lecture_2020_UEA_Short\Computer"
setwd(path_to_folder)

Sigma <- matrix(c(0.1,0,0.23,0,1,0, 0.23,0,1),3,3)


n1 = 7584
#n1 = 758
#n1 = 75

multivariate <- mvrnorm(n = n1, rep(0, 3), Sigma)
multivariate

cor(multivariate[ , 1], multivariate[ ,2])
cor(multivariate[ , 1], multivariate[ ,3])
cor(multivariate[ , 2], multivariate[ ,3])

cov(multivariate[ , 1], multivariate[ ,2])
cov(multivariate[ , 1], multivariate[ ,3])
cov(multivariate[ , 3], multivariate[ ,2])


names(multivariate)[1] <- "x1"
names(multivariate)[2] <- "x2"
#names(multivariate)[3] <- "x3"


summary(multivariate)


write.table(multivariate, "3varsData.RData")
multivarfile = file.path(path_to_folder, "3varsData.RData")
write.table(multivariate, multivarfile)


#load("mvdata.RData")
mvdata <- read.table("3varsData.RData")

cor(mvdata[ , 1], mvdata[ ,2])
cov(mvdata[ , 1], mvdata[ ,2])

library(ggplot2)
qplot (V1, V2, geom="point", data=mvdata)
qplot (V1, V2, geom="point", data=mvdata)+ geom_smooth(method=lm, se=FALSE)
qplot (V1, V3, geom="point", data=mvdata)+ geom_smooth(method=lm, se=FALSE)



sigma_iid <- matrix(c(2,0,0, 0,200,0, 0,0,10),3,3)
eps <- mvrnorm(n = n1, rep(0, 3), sigma_iid)

eps1 = eps[,1]
eps2 = eps[,2]

names(eps)[1] <- "eps1"
names(eps)[2] <- "eps2"
names(eps)[3] <- "eps0"


summary(mvdata)


names(mvdata)[1] <- "x1"
names(mvdata)[2] <- "x2"
names(mvdata)[3] <- "x3"

summary(mvdata)


fivevars = cbind(mvdata, eps[,1], eps[,2])
dim(fivevars)

names(fivevars)[4] <- "eps1"
names(fivevars)[5] <- "eps2"

summary(fivevars)

const = matrix(1, n1, 1)
constant = const*2.2
eps0 <- rnorm(n=n1, 0, 1)

constant = const*2.2
y1 =  constant + 0.38772*fivevars[1] + 1.3221*fivevars[2] + 0.56*fivevars[3] + eps[, 3]

dim(y1) 
dim(transpose(y1))
dim(eps[3])
dim(const)

"sevenvars"
sevenvars = cbind(y1, const, fivevars)
dim(sevenvars)
names(sevenvars)[1]="y1"
names(sevenvars)[2]="const"
summary()


write.table(sevenvars, file="OLSqData.RData")
df.OLSdat <- read.table("OLSqData.RData")
DT <- data.table(df.OLSdat) 

summary(df.OLSdat)
summary(DT)
View(DT)


cor(DT$y1, DT$x1)
cor(DT)
cov(DT)


y1
DT$y1

mod <- lm(formula = DT$y1 ~  DT$x1 + DT$x2 + DT$x3)
summary(mod)




qplot (DT$y1, DT$x2, geom="point", data=mvdata)+ geom_smooth(method=lm, se=FALSE)
qplot (DT$y1, DT$x3, geom="point", data=mvdata)+ geom_smooth(method=lm, se=FALSE)


# Write CSV
fwrite(DT, "OLSqData.csv")
OLSdat_file = file.path(path_to_folder, "OLSqData.csv")
fwrite(DT, OLSdat_file)


X <- cbind(DT$const, DT$x1, DT$x2, DT$x3)
X

y <- cbind(DT$y1)

Xt <- t(X)
Xt


XtX <- Xt%*%X
XtX

invXtX <- inv(XtX)
invXtX

Xty <- Xt%*%y
Xty

beta <- invXtX%*%Xty
beta