home

The following code is provided as is with no warranties express or implied or statutory or whatever, standard disclaimers apply, use at your own risk, yatta, yatta, yatta. If you find any problems, please make a comment here.

Here's a little Scilab code I wrote to calculate a Variance-Covariance Matrix (MathWorld, NIST).

function [y] = varcovar(m)

   // copy source
   x = m;

   // dimensions
   [nrows,ncols] = size(m);


   // scaling factor
   scale = sqrt(1.0 / (nrows-1));


   // subtract mean from each column
   for i=1:ncols

      x(:,i) = (x(:,i) - mean(x(:,i))) * scale;

   end


   // x'x
   y = x' * x;

endfunction



Example using NIST data


-->x = [4.0 2.0 0.60; 4.2 2.1 0.59; 3.9 2.0 0.58;
 4.3 2.1 0.62; 4.1 2.2 0.63]

 x  =

    4.     2.     0.6
    4.2    2.1    0.59
    3.9    2.     0.58
    4.3    2.1    0.62
    4.1    2.2    0.63

-->varcovar(x)

 ans  =

    0.025      0.0075     0.00175
    0.0075     0.007      0.00135
    0.00175    0.00135    0.00043


Note: If you're concerned with issues of stability and performance, see Knuth, et al.