source('../../Data/data_loaders.R')

CXY = load_appendix_lamb_data()
N = 337 # number of samples
m = 6 # number of predictors

# Part (a):
#
# Extract the coefficients of the full linear model:
# As we have a correlation matrix these predictors are centered and normalized (divided by the standard deviation)
#
XTy = as.matrix( CXY[1:m, m+1] )
XTX = as.matrix( CXY[1:m, 1:m] )

beta_hat = solve( XTX, XTy )

# The VIFs:
#
VIFs = diag(solve(XTX))
print('VIFs')
print(VIFs)

# The eigenvalues:
#
ES = eigen(XTX)
print('Eigenvalues')
print(ES$values)

# The correlation between predictors:
#
C = cor(CXY[, seq(1, m)])
print(C) # notice one correlation r_ij > 0.95
mask = (abs(C)>0.95) & (C!=1)
print('correlations above 0.95')
print(unique(C[mask]))