Tips: Compute the mean and the variance of multidimensionnal data

You can easily compute the mean and the variance matrix of multidimensional data.

Assume we are handling this kind of data

// values (b,g,r,ir)
typedef CArrayVector<double, 4> Spectrum;

repeated in space and time. The data are stored in an array

// array of values
typedef CArray<Spectrum> ArraySpectrum;
ArraySpectrum datait;

and we want to compute at each time the (multidimensional) mean of this data set. This can be used using the following code :

// array of mean values
typedef CArrayPoint<Spectrum> PointSpectrum;
PointSpectrum mut(datait.cols());
for (int t= datait.beginCols(); t< datait.endCols(); ++t)
{
  mut[t] = 0.;
  for (int i=datait_.beginRows(); i< datait_.endRows(); ++i)
  { mut[t] += datait_(i,t);}
  mut[t]/= data.sizeRows();
}

The variance matrix (using numerical correction) can be computed using the following code :

// covariances values (b,g,r,ir)
typedef CArraySquare<double, 4> CovSpectrum;
// array of mean values
typedef CArrayPoint<CovSpectrum> PointCov;
PointSpectrum sigmat(datait.cols());
for (int t= datait.beginCols(); t< datait.endCols(); ++t)
{
  CovSpectrum var; var=0.0;
  Spectrum sum = 0.0;
  for (int i=datait_.beginRows(); i< datait_.endRows(); ++i)
  {
    Spectrum dev;
    sum += (dev = datait(i,t) - mut[t]);
    var += dev * dev.transpose();
  }
  sigmat[t] = (var - ((sum * sum.transpose())/datait.sizeCols()) )/datait.sizeCols();
}

stk++ handles transparently the multidimensional nature of the data.

See also

STK::Stat::mean, STK::Stat::varianceWithFixedMean