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 training data is usually not a problem for neural network regression, and not a problem for tree-based regression (random forest, gradient boosting).

The most common way to analyze training data for multicollinearity is to compute the VIF (variance inflation factor) for each column of the data.

VIF is a value between 1.0 and positive infinity (in weird scenarios, a VIF value could be less than 1.0). Briefly, if all column VIF values are less than about 7.0, the data is probably OK in terms of multicollinearity.

if VIF is close to 1.0, the column is not correlated with other cols
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 a 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.

Some time ago, I put together a demo using Python NumPy and the scikit library. This was easy because NumPy and scikit had all the components available.

I decided to entertain myself by refactoring my Python demo to the C# language. This required a bit of effort because I had to implement a no-touch linear regression class, plus several helper functions, such as a program-defined MatWithoutColumn() to mimic the built-in numpy np.delete() function.

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
. . .

The output of my C# VIF demo program is:

Begin variance inflation factor (VIF) demo using C#

Loading synthetic (20) normal data

First three lines:
  -0.1660   0.4406  -0.9998  -0.3953  -0.7065
   0.0776  -0.1616   0.3704  -0.5911   0.7562
  -0.9452   0.3409  -0.1654   0.1174  -0.7192

Begin VIF analysis
col =  0  |  vif = 1.1979
col =  1  |  vif = 1.4590
col =  2  |  vif = 1.2345
col =  3  |  vif = 1.3025
col =  4  |  vif = 1.2120

Loading synthetic (20) highly multicollinear data
(col[2] = 2.0 * col[0] + col[1] + noise)

First three lines:
  -0.1660   0.4406   0.1096  -0.3953  -0.7065
   0.0776  -0.1616  -0.0045  -0.5911   0.7562
  -0.9452   0.3409  -1.5482   0.1174  -0.7192

Begin VIF analysis
col =  0  |  vif = 25101680.1917
col =  1  |  vif = 5769710.7907
col =  2  |  vif = 30816932.6001
col =  3  |  vif = 1.2937
col =  4  |  vif = 1.1879

End 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.

For my internal linear regressor, I used SGD training with some default parameters and auto-early-exit, so that the regressor just works without tuning. While experimenting, I noticed that SGD training has a lot of trouble dealing with multicollinear data — which is the whole point of checking for multicollinear data in the first place.

An interesting exploration.

In machine learning, you don’t want a relationship between two columns in your training data. But in science fiction movies, you absolutely do want a good relationship between the hero and the main actress.

Two huge lost opportunities because of no actor personality multicollinearity.

Demo program. Replace “lt” (less than), “gt”, “lte”, “gte” with Boolean operator symbols. (My lame blog editor chokes on symbols).

using System;
using System.IO;
using System.Collections.Generic;

namespace VarianceInflationFactor
{
  internal class VarianceInflationFactorProgram
  {
    static void Main(string[] args)
    {
      Console.WriteLine("\nBegin variance inflation " +
        "factor (VIF) demo using C# ");

      Console.WriteLine("\nLoading synthetic (20)" +
        " normal data ");
      string file1 = 
        "..\\..\\..\\Data\\synthetic_train_20.txt";
      double[][] train1 = MatLoad(file1, 
        new int[] { 0, 1, 2, 3, 4 }, ',', "#");
      Console.WriteLine("\nFirst three lines: ");
      for (int i = 0; i "lt" 3; ++i)
        VecShow(train1[i], 4, 9);

      Console.WriteLine("\nBegin VIF analysis ");
      for (int j = 0; j "lt" train1[0].Length; ++j)
      {
        double z = VarInfFactor(train1, j);
        Console.WriteLine("col = " +
          j.ToString().PadLeft(2) + 
          "  |  vif = " + z.ToString("F4"));
      }

      Console.WriteLine("\nLoading synthetic (20)" +
        " highly multicollinear data ");
      Console.WriteLine("(col[2] = 2.0 * col[0] +" +
        " col[1] + noise) ");
      string file2 = 
        "..\\..\\..\\Data\\synthetic_train_20_collinear.txt";
      double[][] train2 = MatLoad(file2,
        new int[] { 0, 1, 2, 3, 4 }, ',', "#");
      Console.WriteLine("\nFirst three lines: ");
      for (int i = 0; i "lt" 3; ++i)
        VecShow(train2[i], 4, 9);

      Console.WriteLine("\nBegin VIF analysis ");
      for (int j = 0; j "lt" train2[0].Length; ++j)
      {
        double z = VarInfFactor(train2, j);
        Console.WriteLine("col = " +
          j.ToString().PadLeft(2) +
          "  |  vif = " + z.ToString("F4"));
      }
  
      Console.WriteLine("\nEnd demo ");
      Console.ReadLine();
    } // Main()

    // ------------------------------------------------------

    static double VarInfFactor(double[][] data, int col)
    {
      // predict col in data[][] using other columns
      double[][] X = MatWithoutCol(data, col);
      double[] y = MatGetColumn(data, col);
      LinearRegressor model = new LinearRegressor();
      model.Train(X, y);
      //model.TrainLeftPinv(X, y);
      double r2 = model.R2(X, y);  // usually in (0.0, 1.0)
      double vif = 1.0 / (1.0 - r2);  // r2 could be 1
      return vif;
    }

    // ------------------------------------------------------
    // helpers for Main(): MatLoad(), VecShow()
    // ------------------------------------------------------

    static double[][] MatLoad(string fn, int[] usecols,
      char sep, string comment)
    {
      List"lt"double[]"gt" result = 
        new List"lt"double[]"gt"();
      string line = "";
      FileStream ifs = new FileStream(fn, FileMode.Open);
      StreamReader sr = new StreamReader(ifs);
      while ((line = sr.ReadLine()) != null)
      {
        if (line.StartsWith(comment) == true)
          continue;
        string[] tokens = line.Split(sep);
        List"lt"double"gt" lst = new List"lt"double"gt"();
        for (int j = 0; j "lt" usecols.Length; ++j)
          lst.Add(double.Parse(tokens[usecols[j]]));
        double[] row = lst.ToArray();
        result.Add(row);
      }
      sr.Close(); ifs.Close();
      return result.ToArray();
    }

    static void VecShow(double[] vec, int dec, int wid)
    {
      for (int i = 0; i "lt" vec.Length; ++i)
        Console.Write(vec[i].ToString("F" + dec).
          PadLeft(wid));
      Console.WriteLine("");
    }

    // ------------------------------------------------------
    // helpers for VarInflationFactor()
    // ------------------------------------------------------

    static double[][] MatWithoutCol(double[][] data, int col)
    {
      int nRows = data.Length;
      int nCols = data[0].Length;

      double[][] result = new double[nRows][];
      for (int i = 0; i "lt" nRows; ++i)
        result[i] = new double[nCols-1];

      int k = 0; // into result
      for (int j = 0; j "lt" nCols; ++j)
      {
        if (j == col) continue;
        for (int i = 0; i "lt" nRows; ++i)
          result[i][k] = data[i][j];
        ++k;
      }
      return result;
    }

    static double[] MatGetColumn(double[][] data, int col)
    {
      int nRows = data.Length;
      int nCols = data[0].Length;
      double[] result = new double[nRows];
      for (int i = 0; i "lt" nRows; ++i)
        result[i] = data[i][col];
      return result;
    }

    // ------------------------------------------------------

  } // class Program

  // ========================================================

  public class LinearRegressor
  {
    public double[] weights;
    public double bias;
    private Random rnd;

    public LinearRegressor(int seed = 1)
    {
      this.weights = new double[0]; // keep compiler happy
      this.bias = 0;
      this.rnd = new Random(seed);
    }

    // ------------------------------------------------------

    public double Predict(double[] x)
    {
      double result = 0.0;
      for (int j = 0; j "lt" x.Length; ++j)
        result += x[j] * this.weights[j];
      result += this.bias;
      return result;
    }

    // ------------------------------------------------------

    public int Train(double[][] trainX, double[] trainY,
      double initRate=0.10, int maxEpochs=10000,
      double noChangeTol=0.00001, int consecutiveNoChange=5)
    {
      // exit when dist(old wts, new wts) lt noChangeTol
      // for consecutiveNoChange times
      int n = trainX.Length;
      int dim = trainX[0].Length;
      this.weights = new double[dim];

      // initialize weights and bias
      // technically not necessary
      double low = -0.01; double hi = 0.01;
      for (int i = 0; i "lt" dim; ++i)
        this.weights[i] = (hi - low) *
          this.rnd.NextDouble() + low;
      this.bias = (hi - low) *
          this.rnd.NextDouble() + low;

      int[] indices = new int[n];  // of train data
      for (int i = 0; i "lt" n; ++i)
        indices[i] = i;

      double[] oldWeightsAndB = new double[dim+1];
      double[] newWeightsAndB = new double[dim+1];
      int countNoChange = 0;
      
      for (int epoch = 0; epoch "lt" maxEpochs; ++epoch)
      {
        Shuffle(indices, this.rnd);

        double lrnRate = 
          initRate / Math.Pow((double)(epoch+1), 0.25);
        // Console.WriteLine(lrnRate);

        for (int j = 0; j "lt" dim; ++j)
          oldWeightsAndB[j] = this.weights[j];
        oldWeightsAndB[dim] = this.bias;

        for (int i = 0; i "lt" n; ++i) // each train item
        {
          int idx = indices[i];
          double[] x = trainX[idx];
          double predY = this.Predict(x);
          double actualY = trainY[idx];

          // decay weights before update
          // not needed because regularization not useful
          // for (int j = 0; j "lt" dim; ++j) // each weight
          //   this.weights[j] *= (1 - decay);

          // update weights and bias
          for (int j = 0; j "lt" dim; ++j) // each weight
            this.weights[j] -= lrnRate *
              (predY - actualY) * x[j];
          this.bias -= lrnRate * (predY - actualY) * 1;
        }

        //// display progress 5 times: very noisy
        //if (epoch % (int)(maxEpochs / 5) == 0) // progress
        //{
        //  double r2 = this.R2(trainX, trainY);
        //  string s = "";
        //  s += "epoch = " + epoch.ToString().PadLeft(5);
        //  s += "  R2 = " + r2.ToString("F4").PadLeft(8);
        //  Console.WriteLine(s);
        //}

        // check for early-exit after each epoch
        for (int j = 0; j "lt" dim; ++j)
          newWeightsAndB[j] = this.weights[j];
        newWeightsAndB[dim] = this.bias;
        double d = EuclideanDist(oldWeightsAndB,
          newWeightsAndB);

        if (d "lt" noChangeTol)
        {
          ++countNoChange;
          if (countNoChange == consecutiveNoChange)
            return epoch;
        }
        else
         countNoChange = 0; // reset

      } // epoch
      return maxEpochs;

    } // Train

    // ------------------------------------------------------

    private static double EuclideanDist(double[] v1,
      double[] v2)
    {
      int n = v1.Length;
      double sum = 0.0;
      for (int i = 0; i "lt" n; ++i)
        sum += (v1[i] - v2[i]) * (v1[i] - v2[i]);
      return Math.Sqrt(sum);
    }

    // ------------------------------------------------------

    public double R2(double[][] dataX, double[] dataY)
    {
      int n = dataX.Length;
      double sum = 0.0;

      for (int i = 0; i "lt" n; ++i)
        sum += dataY[i];
      double meanY = sum / n;

      double ssRes = 0.0; // sum squared residuals
      double ssTot = 0.0; // sum squared total
      for (int i = 0; i "lt" n; ++i)
      {
        double predY = this.Predict(dataX[i]);
        ssRes +=
          (dataY[i] - predY) * (dataY[i] - predY);
        ssTot +=
          (dataY[i] - meanY) * (dataY[i] - meanY);
      }
      return 1.0 - (ssRes / ssTot);
    }

    // ------------------------------------------------------

    private static void Shuffle(int[] indices, Random rnd)
    {
      int n = indices.Length;
      for (int i = 0; i "lt" n; ++i) // one extra pass
      {
        int ri = rnd.Next(i, n);
        int tmp = indices[i];
        indices[i] = indices[ri];
        indices[ri] = tmp;
      }
    }

  } // class LinearRegressor

} // ns

First, normal, no multicollinearity 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, highly 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