Regression : Variance Inflation Factor

[sourish-cmi]

Is your feature request related to a problem? Please describe.
Add Variance Inflation Factor to Frequentist Regression.

Describe the solution you'd like

using DataFrames, GLM, Statistics, LinearAlgebra, RDatasets 
airquality = rename(dataset("datasets", "airquality"), "Solar.R" => "Solar_R")
#
y =  airquality.Wind
X1 = airquality.Temp    
X2 = airquality.Solar_R  
X3 = airquality.Ozone 
d = DataFrame(X1=X1, X2=X2, X3=X3, y=y)
d = d[completecases(d), :]
X = [one.(d.X1) d.X1 d.X2 d.X3]  .|> Float64
y = d.y .|> Float64
#
β_LS = X \ y
vifm = diag(inv(cor(X[:,2:end])))
#
m_y  = lm(@formula(y ~ 1 + X1 + X2 + X3), d)
m_X1 = lm(@formula(X1 ~ 1 + X2 + X3), d)
m_X2 = lm(@formula(X2 ~ 1 + X1 + X3), d)
m_X3 = lm(@formula(X3 ~ 1 + X1 + X2), d)

vif_X1 = 1.0/(1.0-r2(m_X1))
vif_X2 = 1.0/(1.0-r2(m_X2))
vif_X3 = 1.0/(1.0-r2(m_X3))
#
vifm ≈ [vif_X1; vif_X2; vif_X3;]
#
vif_GLM(mod) = diag(inv(cor(mod.model.pp.X[:,2:end])))
vif_GLM(m_y) ≈ vifm ≈ [vif_X1; vif_X2; vif_X3;]

Detailed discussion can be found here

The function should be like below:

vif_hat = vif(@formula(y ~ 1 + X1 + X2 + X3), d)

vif_hat returns vifm ≈ [vif_X1; vif_X2; vif_X3;]

Describe alternatives you've considered
If it exists in StatsAPI or any other native Julia package, integrate it

Additional context
nothing as such