Distance calculation between multi-dimension gaussian distributions (diagonal co-variance). 多次元正規分布（対角共分散）間の距離計算方法

CalcGaussianDistance.cpp

#include <cmath>

// ----------------------------------------------------------------------------
// Kullback-Leiblier Divergence (KLD-divergence, asymmetric)
// ----------------------------------------------------------------------------
float calcKld(
    const int    nDim,
	const float* mean1,
	const float* var1,
	const float* mean2,
	const float* var2
)
{
	float	sum = 0.0;
	float	term1 = 0.0;
	float	term2 = 0.0;
	float	term3 = 0.0;

	for( int d = 0; d < nDim; ++d )
	{
		term1 = var1[ d ] / var2[ d ];
		term2 = ( mean2[ d ] - mean1[ d ] ) * ( mean2[ d ] - mean1[ d ] ) / var2[ d ];
		term3 = log( var1[ d ] / var2[ d ] );
		sum += ( term1 + term2 - term3 - 1 );
	}

	return sum * 0.5f;
}

// ----------------------------------------------------------------------------
// KLD Average
// ----------------------------------------------------------------------------
float calcKldAvg(
	const int	nDim,
	const float* mean1,
	const float* var1,
	const float* mean2,
	const float* var2
)
{
	float	sum = 0.0;
	float	meanDiffSq = 0.0;
	float	term1 = 0.0;
	float	term2 = 0.0;

	for( int d = 0; d < nDim; ++d )
	{
		meanDiffSq	= ( mean2[ d ] - mean1[ d ] ) * ( mean2[ d ] - mean1[ d ] );
		term1		= ( var1[ d ] + meanDiffSq ) / var2[ d ];
		term2		= ( var2[ d ] + meanDiffSq ) / var1[ d ];
		sum			+= ( term1 + term2 - 2 );
	}

	return sum * 0.25f;
}


// ----------------------------------------------------------------------------
// Bhatacharyya Distance (asymmetric)
// ----------------------------------------------------------------------------
float calcBhattacharyyaDistance(
	const int	nDim,
	const float* mean1,
	const float* var1,
	const float* mean2,
	const float* var2
)
{
	float	sum = 0.0;
	float	term1 = 0.0;
	float	term2 = 0.0;

	for( int d = 0; d < nDim; ++d )
	{
		term1 = ( mean1[ d ] - mean2[ d ] ) * ( mean1[ d ] - mean2[ d ] ) / ( var1[ d ] + var2[ d ] );
		term2 = log( ( var1[ d ] + var2[ d ] ) * ( var1[ d ] + var2[ d ] ) / ( 4 * var1[ d ] * var2[ d ] ) );
		sum += term1 + term2;
	}

	return sum * 0.25f;
}


// ----------------------------------------------------------------------------
// Euclidean Distance (symmetric)
// ----------------------------------------------------------------------------
float calcEuclideanDistance(
	const int	nDim,
	const float* mean1,
	const float* var1,
	const float* mean2,
	const float* var2
)
{
	float	sum = 0.0;
	float	diff = 0.0;

	for( int d = 0; d < nDim; ++d )
	{
		diff = mean1[ d ] - mean2[ d ];
		sum += diff * diff;
	}

	return sqrt( sum );
}