In machine learning, if training data is multicollinear, the interpretability of the resulting model will likely be poor (the details are tricky and out of the scope of this post).
Multicollinear data has two or more rows that are highly correlated mathematically. For example, if a set of training data has a column (predictor variable) house-size (in square feet), and another column house-price, those two columns are likely highly correlated with each other.
Note: Multicollinear data is usually not a problem for neural network regression, and not a problem for tree-based regression (random forest, gradient boosting).
I have implemented VIF analysis using from-scratch C# and from-scratch Python. One morning before work, I figured I’d implement VIF analysis from-scratch JavaScript.
VIF is a value between 1.0 and positive infinity (actually, in weird scenarios, a VIF value could be less than one). Briefly, if all column VIF values are less than about 7.0, the data is probably OK.
if VIF is close to 1.0, the column is not correlated.
if VIF between 1.0 and 5.0, column is mildly correlated
if VIF between 5.0 and 10.0, column is highly correlated
if VIF greater than 10.0, column is extremely correlated
To compute the VIF for a specified column of training data, you use the specified column as the dependent y variable, and use the remaining columns as the independent predictor variables, and compute a linear regression model, and then compute the R2 (coefficient of determination) for the model. The VIF value for the column is 1.0 / (1.0 – R2).
Suppose that you have a set of training data X predictor values, and you use some column c as the dependent y variable, and all the other columns as predictors for c. After training the linear regression model, you compute R2 and it is 0.90 — which means column c is predicted very well by the other columns. The VIF value for column c is 1.0 / (1.0 – R2) = 1.0 / 0.10 = 10.0 which is large which is bad because column c is a linear combination of the other columns — the data is somewhat multicollinear.
Now, with the same setup, suppose R2 is 0.20 — which means column c cannot be predicted well by the other columns. The VIF value is 1.0 / (1.0 – 0.20) = 1.0 / 0.8 = 1.25 which is a small value, which is good, because column c is not a linear combination of the other columns, and therefore the data is not multicollinear.
I created two datasets. The first data set has five columns of predictors, followed by a column of target y values. The data is “normal” in the sense that there’s no multicollinearity. There are 20 items. It looks like:
-0.1660, 0.4406, -0.9998, -0.3953, -0.7065, 0.4840 0.0776, -0.1616, 0.3704, -0.5911, 0.7562, 0.1568 -0.9452, 0.3409, -0.1654, 0.1174, -0.7192, 0.8054 . . .
The second dataset is highly multicollinear, where the third column is 2 times the first column, plus the second column, plus a small random value between 0.000 and 0.001. It looks like:
-0.1660, 0.4406, 0.1096, -0.3953, -0.7065, 0.4840 0.0776, -0.1616, -0.0045, -0.5911, 0.7562, 0.1568 -0.9452, 0.3409, -1.5482, 0.1174, -0.7192, 0.8054 . . .
After about 45 minutes of coding, I had a JavaScript program up and running. The output of the demo program is:
Begin variance inflation factor (VIF) demo using JavaScript Loading synthetic (20) normal non-collinear dataset from file First two items: -0.1660 0.4406 -0.9998 -0.3953 -0.7065 0.0776 -0.1616 0.3704 -0.5911 0.7562 Begin VIF analysis col = 0 | vif = 1.1980 col = 1 | vif = 1.4591 col = 2 | vif = 1.2345 col = 3 | vif = 1.3025 col = 4 | vif = 1.2120 Loading synthetic (20) highly collinear dataset from file (col[2] = 2.0 * col[0] + col[1] + rnd) First two items: -0.1660 0.4406 0.1096 -0.3953 -0.7065 0.0776 -0.1616 -0.0045 -0.5911 0.7562 Begin VIF analysis col = 0 | vif = 25546262.9389 col = 1 | vif = 6023299.3846 col = 2 | vif = 30951889.0306 col = 3 | vif = 1.2937 col = 4 | vif = 1.2117 End VIF demo
As expected, the first dataset didn’t have any bad VIF values, but the VIF values for the second dataset show that columns [0], [1], [2] are highly correlated.
Behind the scenes, I train each linear regression model using a closed-form algorithm — relaxed Moore-Penrose pseudo-inverse via QR decomposition (Householder version). This could fail for very large datasets. In such cases, you can use stochastic gradient descent (SGD), which can handle arbitrarily large datasets.
No big moral to this blog post. Just an interesting exploration.
I’m a big fan machine learning, and old science fiction movies. In machine learning, highly correlated data columns are a bad thing. In old science fiction movies, highly correlated couples are a good thing.
Left: In “Earth vs. the Flying Saucers” (1956), scientist Dr. Russell Marvin and his new wife Carol, have good chemistry, even when they’re being stalked by a flying saucer, visible in the rear window of their car.
Right: In “Beginning of the End” (1957), scientist Dr. Ed Wainwright and romantic interest, news reporter Audrey Aimes, have good chemistry, even when they’re about to encounter huge, man-eating grasshoppers.
Demo program. Replace “lt” (lss than), “gt”, “lte”, “gte” with Boolean operator symbols (my blog editor chokes on symbols).
// variance_inflation_factor.js
// compute VIF: are columns of a dataset collinear?
let FS = require("fs") // for loadTxt()
// ==========================================================
function varInfFactor(data, col)
{
// treat col as dependent variable, other columns as preds
let X = matDeleteCol(data, col);
let y = matGetColumn(data, col)
let model = new LinearRegressor();
model.train(X, y);
let r2 = model.R2(X, y);
let vif = 1.0 / (1.0 - r2);
return vif;
}
// helpers for varInfFactor()
function matDeleteCol(data, col)
{
let nRows = data.length;
let nCols = data[0].length;
let result = matMake(nRows, nCols-1, 0.0);
let k = 0; // into result
for (let j = 0; j "lt" nCols; ++j) {
if (j == col) continue;
for (let i = 0; i "lt" nRows; ++i)
result[i][k] = data[i][j];
++k;
}
return result;
}
function matGetColumn(data, col)
{
let nRows = data.length;
let nCols = data[0].length;
let result = vecMake(nRows, 0.0);
for (let i = 0; i "lt" nRows; ++i)
result[i] = data[i][col];
return result;
}
// minor helpers for the helpers
function matMake(nRows, nCols, val)
{
let result = [];
for (let i = 0; i "lt" nRows; ++i) {
result[i] = [];
for (let j = 0; j "lt" nCols; ++j) {
result[i][j] = val;
}
}
return result;
}
function vecMake(n, val)
{
let result = [];
for (let i = 0; i "lt" n; ++i) {
result[i] = val;
}
return result;
}
// ==========================================================
class LinearRegressor
{
// LR using closed-form MP pseudo-inverse training
constructor(seed)
{
this.weights; // allocated in train()
this.bias = 0.0; // supplied in train()
}
// --------------------------------------------------------
// methods: train(), predict(), R2()
// --------------------------------------------------------
predict(x)
{
let sum = 0.0;
for (let i = 0; i "lt" x.length; ++i) {
sum += x[i] * this.weights[i];
}
sum += this.bias;
return sum;
}
// --------------------------------------------------------
train(trainX, trainY)
{
// w = pinv(DX) * y
let dim = trainX[0].length; // number predictors
this.weights = this.vecMake(dim, 0.0); // allocate wts
let DX = this.matToDesign(trainX); // design matrix
let Xpinv = this.matPinv(DX); // pinv of design X
let biasAndWts = this.matVecProd(Xpinv, trainY);
this.bias = biasAndWts[0];
for (let i = 1; i "lt" biasAndWts.length; ++i)
this.weights[i-1] = biasAndWts[i];
return;
}
// --------------------------------------------------------
R2(dataX, dataY)
{
let n = dataX.length;
let sum = 0.0;
for (let i = 0; i "lt" n; ++i)
sum += dataY[i];
let meanY = sum / n;
let ssRes = 0.0; // sum squared residuals
let ssTot = 0.0; // sum squared total
for (let i = 0; i "lt" n; ++i) {
let predY = this.predict(dataX[i]);
ssRes += (dataY[i] - predY) * (dataY[i] - predY);
ssTot += (dataY[i] - meanY) * (dataY[i] - meanY);
}
return 1.0 - (ssRes / ssTot);
}
// --------------------------------------------------------
// primary helpers for train()
// matVecProd, matPinv, matToDesign
// --------------------------------------------------------
matVecProd(M, v)
{
// return a regular 1D vector
let nRows = M.length;
let nCols = M[0].length;
let n = v.length; // assume nCols = n
let result = this.vecMake(nRows, 0.0);
for (let i = 0; i "lt" nRows; ++i)
for (let k = 0; k "lt" nCols; ++k)
result[i] += M[i][k] * v[k];
return result;
}
// --------------------------------------------------------
matToDesign(M)
{
// add a leading column of 1.0s to M
let nRows = M.length;
let nCols = M[0].length;
let result = this.matMake(nRows, nCols+1, 1.0);
for (let i = 0; i "lt" nRows; ++i) {
for (let j = 1; j "lt" nCols+1; ++j) { // note 1s
result[i][j] = M[i][j-1];
}
}
return result;
}
// --------------------------------------------------------
matPinv(M) // Moore-Penrose using QR-Householder decomp
{
// A = Q*R, pinv(A) = inv(R) * trans(Q)
let m = M.length; let n = M[0].length;
if (m "lt" n)
cosole.log("ERROR. works only m "gte" m");
let QR = this.matDecomposeQR(M, true); // reduced
let Rinv = this.matInvUpperTri(QR[1]);
let Qinv = this.matTranspose(QR[0]);
let result = this.matProduct(Rinv, Qinv);
return result;
}
// --------------------------------------------------------
// helpers for matPinv
// matDecomposeQR, matInvUpperTri, matTranspose, matProduct
// --------------------------------------------------------
matDecomposeQR(M, reduced)
{
// QR decomposition, Householder algorithm.
// result[0] = Q, result[1] = R
let m = M.length;
let n = M[0].length;
if (m "lt" n)
console.log("No rows less than cols");
let Q = this.matIdentity(m); // working Q
let R = this.matCopy(M); // working R
let end = 0;
if (m == n) end = n - 1;
else end = n;
for (let i = 0; i "lt" end; ++i) {
let H = this.matIdentity(m);
let a = this.vecMake(m-i, 0.0);
let k = 0;
for (let ii = i; ii "lt" m; ++ii) // corr
a[k++] = R[ii][i];
let normA = this.vecNorm(a);
if (a[0] "lt" 0.0 && normA "gt" 0.0) // corr
normA = -normA;
else if (a[0] "gt" 0.0 && normA "lt" 0.0)
normA = -normA;
let v = this.vecMake(a.length, 0.0);
for (let j = 0; j "lt" v.length; ++j)
v[j] = a[j] / (a[0] + normA);
v[0] = 1.0;
// Householder algorithm
let h = this.matIdentity(a.length);
let vvDot = this.vecDot(v, v);
let A = this.vecToMat(v, v.length, 1);
let B = this.vecToMat(v, 1, v.length);
let AB = this.matProduct(A, B);
for (let ii = 0; ii "lt" h.length; ++ii)
for (let jj = 0; jj "lt" h[0].length; ++jj)
h[ii][jj] -= (2.0 / vvDot) * AB[ii][jj];
// copy h[][] into lower right corner of H[][]
let d = m - h.length; // corr
for (let ii = 0; ii "lt" h.length; ++ii)
for (let jj = 0; jj "lt" h[0].length; ++jj)
H[ii + d][jj + d] = h[ii][jj];
Q = this.matProduct(Q, H);
R = this.matProduct(H, R);
} // i
let result = [];
if (reduced == false) {
result[0] = Q;
result[1] = R;
return result;
}
else if (reduced == true) {
let qRows = Q.length; let qCols = Q[0].length;
let rRows = R.length; let rCols = R[0].length;
// assumes m "gte" n !!
// square-up R
let dim = Math.min(rRows, rCols);
let Rsquared = this.matMake(dim, dim, 0.0);
for (let i = 0; i "lt" dim; ++i)
for (let j = 0; j "lt" dim; ++j)
Rsquared[i][j] = R[i][j];
// Q needs same number columns as R
// so that inv(R) * trans(Q) works
let Qtrimmed = this.matMake(qRows, dim, 0.0);
for (let i = 0; i "lt" qRows; ++i)
for (let j = 0; j "lt" dim; ++j)
Qtrimmed[i][j] = Q[i][j];
result[0] = Qtrimmed;
result[1] = Rsquared;
return result;
}
} // matDecomposeQR()
// --------------------------------------------------------
matInvUpperTri(U)
{
let n = U.length; // U must be square matrix
let result = this.matIdentity(n);
for (let k = 0; k "lt" n; ++k) {
for (let j = 0; j "lt" n; ++j) {
for (let i = 0; i "lt" k; ++i) {
result[j][k] -= result[j][i] * U[i][k];
}
result[j][k] /= (U[k][k] + 1.0e-8); // avoid 0
}
}
return result;
}
// --------------------------------------------------------
matTranspose(M)
{
let nRows = M.length;
let nCols = M[0].length;
let result = this.matMake(nCols, nRows, 0.0);
for (let i = 0; i "lt" nRows; ++i)
for (let j = 0; j "lt" nCols; ++j)
result[j][i] = M[i][j];
return result;
}
// --------------------------------------------------------
matProduct(matA, matB)
{
let aRows = matA.length;
let aCols = matA[0].length;
let bRows = matB.length;
let bCols = matB[0].length;
let result = this.matMake(aRows, bCols, 0.0);
for (let i = 0; i "lt" aRows; ++i)
for (let j = 0; j "lt" bCols; ++j)
for (let k = 0; k "lt" aCols; ++k)
result[i][j] += matA[i][k] * matB[k][j];
return result;
}
// --------------------------------------------------------
// misc helper methods/functions:
// matMake, matCopy, matIdentity, vecMake,
// vecToMat, vecNorm, vecDot
//
// all of these could be declared 'static'
// --------------------------------------------------------
matMake(nRows, nCols, val)
{
let result = [];
for (let i = 0; i "lt" nRows; ++i) {
result[i] = [];
for (let j = 0; j "lt" nCols; ++j) {
result[i][j] = val;
}
}
return result;
}
// --------------------------------------------------------
matCopy(M)
{
let nRows = M.length;
let nCols = M[0].length;
let result = this.matMake(nRows, nCols, 0.0);
for (let i = 0; i "lt" nRows; ++i)
for (let j = 0; j "lt" nCols; ++j)
result[i][j] = M[i][j];
return result;
}
// --------------------------------------------------------
matIdentity(n)
{
let result = this.matMake(n, n, 0.0);
for (let i = 0; i "lt" n; ++i)
result[i][i] = 1.0;
return result;
}
// --------------------------------------------------------
vecMake(n, val)
{
let result = [];
for (let i = 0; i "lt" n; ++i) {
result[i] = val;
}
return result;
}
// --------------------------------------------------------
vecToMat(vec, nRows, nCols)
{
let result = this.matMake(nRows, nCols, 0.0);
let k = 0;
for (let i = 0; i "lt" nRows; ++i)
for (let j = 0; j "lt" nCols; ++j)
result[i][j] = vec[k++];
return result;
}
// --------------------------------------------------------
vecNorm(vec)
{
let sum = 0.0;
let n = vec.length;
for (let i = 0; i "lt" n; ++i)
sum += vec[i] * vec[i];
return Math.sqrt(sum);
}
// --------------------------------------------------------
vecDot(v1, v2)
{
let n = v1.length; // assume len(v1) == len(v2)
let sum = 0.0;
for (let i = 0; i "lt" n; ++i)
sum += v1[i] * v2[i];
return sum;
}
// --------------------------------------------------------
// --------------------------------------------------------
// member methods for debugging: matShow, vecShow
// --------------------------------------------------------
matShow(M, dec, wid) // for debugging
{
let small = 1.0 / Math.pow(10, dec);
let nr = M.length;
let nc = M[0].length;
for (let i = 0; i "lt" nr; ++i) {
for (let j = 0; j "lt" nc; ++j) {
let x = M[i][j];
if (Math.abs(x) "lt" small) x = 0.0;
let xx = x.toFixed(dec);
let s = xx.toString().padStart(wid, ' ');
process.stdout.write(s);
process.stdout.write(" ");
}
process.stdout.write("\n");
}
}
// --------------------------------------------------------
vecShow(vec, dec, wid) // for debugging
{
let small = 1.0 / Math.pow(10, dec);
for (let i = 0; i "lt" vec.length; ++i) {
let x = vec[i];
if (Math.abs(x) "lt" small) x = 0.0 // avoid -0.00
let xx = x.toFixed(dec);
let s = xx.toString().padStart(wid, ' ');
process.stdout.write(s);
process.stdout.write(" ");
}
process.stdout.write("\n");
}
} // end class LinearRegressor
// ==========================================================
// ----------------------------------------------------------
// helper functions for main()
// matShow, vecShow, matLoad, matToVec
// ----------------------------------------------------------
function vecShow(vec, dec, wid, nl)
{
let small = 1.0 / Math.pow(10, dec);
for (let i = 0; i "lt" vec.length; ++i) {
let x = vec[i];
if (Math.abs(x) "lt" small) x = 0.0 // avoid -0.00
let xx = x.toFixed(dec);
let s = xx.toString().padStart(wid, ' ');
process.stdout.write(s);
process.stdout.write(" ");
}
if (nl == true)
process.stdout.write("\n");
}
// ----------------------------------------------------------
function matShow(A, dec, wid)
{
let small = 1.0 / Math.pow(10, dec);
let nr = A.length;
let nc = A[0].length;
for (let i = 0; i "lt" nr; ++i) {
for (let j = 0; j "lt" nc; ++j) {
let x = A[i][j];
if (Math.abs(x) "lt" small) x = 0.0;
let xx = x.toFixed(dec);
let s = xx.toString().padStart(wid, ' ');
process.stdout.write(s);
process.stdout.write(" ");
}
process.stdout.write("\n");
}
}
// ----------------------------------------------------------
function matToVec(M)
{
let nr = M.length;
let nc = M[0].length;
let result = [];
for (let i = 0; i "lt" nr*nc; ++i) { // vecMake(r*c, 0.0);
result[i] = 0.0;
}
let k = 0;
for (let i = 0; i "lt" nr; ++i) {
for (let j = 0; j "lt" nc; ++j) {
result[k++] = M[i][j];
}
}
return result;
}
// ----------------------------------------------------------
function matLoad(fn, delimit, usecols, comment)
{
// efficient but mildly complicated
let all = FS.readFileSync(fn, "utf8"); // giant string
all = all.trim(); // strip final crlf in file
let lines = all.split("\n"); // array of lines
// count number non-comment lines
let nRows = 0;
for (let i = 0; i "lt" lines.length; ++i) {
if (!lines[i].startsWith(comment))
++nRows;
}
let nCols = usecols.length;
// let result = matMake(nRows, nCols, 0.0);
let result = [];
for (let i = 0; i "lt" nRows; ++i) {
result[i] = [];
for (let j = 0; j "lt" nCols; ++j) {
result[i][j] = 0.0;
}
}
let r = 0; // into lines
let i = 0; // into result[][]
while (r "lt" lines.length) {
if (lines[r].startsWith(comment)) {
++r; // next row
}
else {
let tokens = lines[r].split(delimit);
for (let j = 0; j "lt" nCols; ++j) {
result[i][j] = parseFloat(tokens[usecols[j]]);
}
++r;
++i;
}
}
return result;
}
// ==========================================================
function main()
{
console.log("\nBegin variance inflation factor (VIF) " +
"demo using JavaScript ");
// non-collinear (good) dataset
console.log("\nLoading synthetic (20) normal " +
"non-collinear dataset from file ");
let file1 = ".\\Data\\synthetic_train_20.txt";
let train1 = matLoad(file1, ",", [0,1,2,3,4], "#");
console.log("\nFirst two items: ");
for (let i = 0; i "lt" 2; ++i)
vecShow(train1[i], 4, 9, true);
console.log("\nBegin VIF analysis ");
for (let j = 0; j "lt" train1[0].length; ++j) {
let z = varInfFactor(train1, j);
console.log("col = " + j.toString().padStart(2) +
" | vif = " + z.toFixed(4).toString() );
}
// collinear (bad) dataset
console.log("\nLoading synthetic (20) highly collinear " +
"dataset from file ");
console.log("(col[2] = 2.0 * col[0] + col[1] + rnd) " );
let file2 = ".\\Data\\synthetic_train_20_collinear.txt";
let train2 = matLoad(file2, ",", [0,1,2,3,4], "#");
console.log("\nFirst two items: ");
for (let i = 0; i "lt" 2; ++i)
vecShow(train2[i], 4, 9, true);
console.log("\nBegin VIF analysis ");
for (let j = 0; j "lt" train2[0].length; ++j) {
let z = varInfFactor(train2, j);
console.log("col = " + j.toString().padStart(2) +
" | vif = " + z.toFixed(4).toString());
}
console.log("\nEnd VIF demo");
}
main();
First, normal, dataset:
# synthetic_train_20.txt # -0.1660, 0.4406, -0.9998, -0.3953, -0.7065, 0.4840 0.0776, -0.1616, 0.3704, -0.5911, 0.7562, 0.1568 -0.9452, 0.3409, -0.1654, 0.1174, -0.7192, 0.8054 0.9365, -0.3732, 0.3846, 0.7528, 0.7892, 0.1345 -0.8299, -0.9219, -0.6603, 0.7563, -0.8033, 0.7955 0.0663, 0.3838, -0.3690, 0.3730, 0.6693, 0.3206 -0.9634, 0.5003, 0.9777, 0.4963, -0.4391, 0.7377 -0.1042, 0.8172, -0.4128, -0.4244, -0.7399, 0.4801 -0.9613, 0.3577, -0.5767, -0.4689, -0.0169, 0.6861 -0.7065, 0.1786, 0.3995, -0.7953, -0.1719, 0.5569 0.3888, -0.1716, -0.9001, 0.0718, 0.3276, 0.2500 0.1731, 0.8068, -0.7251, -0.7214, 0.6148, 0.3297 -0.2046, -0.6693, 0.8550, -0.3045, 0.5016, 0.2129 0.2473, 0.5019, -0.3022, -0.4601, 0.7918, 0.2613 -0.1438, 0.9297, 0.3269, 0.2434, -0.7705, 0.5171 0.1568, -0.1837, -0.5259, 0.8068, 0.1474, 0.3307 -0.9943, 0.2343, -0.3467, 0.0541, 0.7719, 0.5581 0.2467, -0.9684, 0.8589, 0.3818, 0.9946, 0.1092 -0.6553, -0.7257, 0.8652, 0.3936, -0.8680, 0.7018 0.8460, 0.4230, -0.7515, -0.9602, -0.9476, 0.1996
Second, multicollinear, dataset:
# synthetic_train_20_collinear.txt # col [2] = 2*[0] + [1] + rand(0.001) # -0.1660, 0.4406, 0.1096, -0.3953, -0.7065, 0.4840 0.0776, -0.1616, -0.0045, -0.5911, 0.7562, 0.1568 -0.9452, 0.3409, -1.5482, 0.1174, -0.7192, 0.8054 0.9365, -0.3732, 1.5016, 0.7528, 0.7892, 0.1345 -0.8299, -0.9219, -2.5800, 0.7563, -0.8033, 0.7955 0.0663, 0.3838, 0.5179, 0.3730, 0.6693, 0.3206 -0.9634, 0.5003, -1.4245, 0.4963, -0.4391, 0.7377 -0.1042, 0.8172, 0.6100, -0.4244, -0.7399, 0.4801 -0.9613, 0.3577, -1.5636, -0.4689, -0.0169, 0.6861 -0.7065, 0.1786, -1.2325, -0.7953, -0.1719, 0.5569 0.3888, -0.1716, 0.6073, 0.0718, 0.3276, 0.2500 0.1731, 0.8068, 1.1544, -0.7214, 0.6148, 0.3297 -0.2046, -0.6693, -1.0770, -0.3045, 0.5016, 0.2129 0.2473, 0.5019, 0.9980, -0.4601, 0.7918, 0.2613 -0.1438, 0.9297, 0.6435, 0.2434, -0.7705, 0.5171 0.1568, -0.1837, 0.1313, 0.8068, 0.1474, 0.3307 -0.9943, 0.2343, -1.7528, 0.0541, 0.7719, 0.5581 0.2467, -0.9684, -0.4732, 0.3818, 0.9946, 0.1092 -0.6553, -0.7257, -2.0345, 0.3936, -0.8680, 0.7018 0.8460, 0.4230, 2.1166, -0.9602, -0.9476, 0.1996
.NET Test Automation Recipes
Software Testing
SciPy Programming Succinctly
Keras Succinctly
R Programming
Visual Studio Live
Microsoft MLADS Conference
DevIntersection Conference
Machine Learning Week
Ai4 Conference
G2E Conference
iSC West Conference
You must be logged in to post a comment.