#!/usr/bin/perl
use warnings;
use strict;
use PDL;

# ================================
# covariance: 
#
#   $Sigma = covariance( $X )
#
#   computes the Sample Covariance Matrix of
#   a sample X1...Xn of p-dimensional vectors
# ================================
sub covariance {
    my ( $X ) = @_;
    
    my $Diff = $X - average( $X->xchg(0,1) );
    
    my $Sigma = ( 1 / ( $X->getdim(1) - 1 ) )
                * transpose( $Diff ) x $Diff;
    
    return $Sigma;
}