From 376bb903247016ef0d4eb293c35bd38a312b4c1e Mon Sep 17 00:00:00 2001 From: Moritz Sallermann Date: Sun, 24 Mar 2024 20:20:33 +0000 Subject: [PATCH] Test: Gaussian copula Fixed issues with bivariate_gaussian_copula, due to previous renaming of cdf. Also put it into the same unit test and plotted the histograms. Co-authored-by: Amrita Goswami --- include/util/math.hpp | 15 ++++++------ test/test_probability_distributions.cpp | 31 ++++++------------------- 2 files changed, 14 insertions(+), 32 deletions(-) diff --git a/include/util/math.hpp b/include/util/math.hpp index 4d49810..493e6ab 100644 --- a/include/util/math.hpp +++ b/include/util/math.hpp @@ -243,11 +243,11 @@ class bivariate_gaussian_copula { private: ScalarT covariance; - bivariate_normal_distribution bivariate_normal_dist{}; + bivariate_normal_distribution biv_normal_dist{}; // std::normal_distribution normal_dist{}; // Cumulative probability function for gaussian with mean 0 and variance 1 - ScalarT cum_gauss( ScalarT x ) + ScalarT cdf_gauss( ScalarT x ) { return 0.5 * ( 1 + std::erf( ( x ) / std::sqrt( 2.0 ) ) ); } @@ -258,9 +258,9 @@ class bivariate_gaussian_copula public: bivariate_gaussian_copula( ScalarT covariance, dist1T dist1, dist2T dist2 ) : covariance( covariance ), + biv_normal_dist( bivariate_normal_distribution( covariance ) ), dist1( dist1 ), - dist2( dist2 ), - bivariate_normal_dist( bivariate_normal_dist( covariance ) ) + dist2( dist2 ) { } @@ -268,12 +268,11 @@ class bivariate_gaussian_copula std::array operator()( Generator & gen ) { // 1. Draw from bivariate gaussian - auto z = bivariate_normal_dist( gen ); + auto z = biv_normal_dist( gen ); // 2. Transform marginals to unit interval - std::array z_unit = { cum_gauss( z[0] ), cum_gauss( z[1] ) }; + std::array z_unit = { cdf_gauss( z[0] ), cdf_gauss( z[1] ) }; // 3. Apply inverse transform sampling - std::array res - = { dist1.inverse_cumulative_probability( z_unit[0] ), dist2.inverse_cumulative_probability( z_unit[1] ) }; + std::array res = { dist1.inverse_cdf( z_unit[0] ), dist2.inverse_cdf( z_unit[1] ) }; return res; } }; diff --git a/test/test_probability_distributions.cpp b/test/test_probability_distributions.cpp index 930e4d6..e72a3df 100644 --- a/test/test_probability_distributions.cpp +++ b/test/test_probability_distributions.cpp @@ -46,27 +46,10 @@ TEST_CASE( "Test the probability distributions", "[prob]" ) write_results_to_file( 10000, Seldon::bivariate_normal_distribution( 0.5 ), "bivariate_normal.txt" ); } -// TEST_CASE( "Test reading in the agents from a file", "[io_agents]" ) -// { -// using namespace Seldon; -// using namespace Catch::Matchers; - -// auto proj_root_path = fs::current_path(); -// auto network_file = proj_root_path / fs::path( "test/res/opinions.txt" ); - -// auto agents = Seldon::agents_from_file( network_file ); - -// std::vector opinions_expected = { 2.1127107987061544, 0.8088982488089491, -0.8802809369462433 }; -// std::vector activities_expected = { 0.044554683389757696, 0.015813166022685163, 0.015863953902810535 }; -// std::vector reluctances_expected = { 1.0, 1.0, 2.3 }; - -// REQUIRE( agents.size() == 3 ); - -// for( size_t i = 0; i < agents.size(); i++ ) -// { -// fmt::print( "{}", i ); -// REQUIRE_THAT( agents[i].data.opinion, Catch::Matchers::WithinAbs( opinions_expected[i], 1e-16 ) ); -// REQUIRE_THAT( agents[i].data.activity, Catch::Matchers::WithinAbs( activities_expected[i], 1e-16 ) ); -// REQUIRE_THAT( agents[i].data.reluctance, Catch::Matchers::WithinAbs( reluctances_expected[i], 1e-16 ) ); -// } -// } \ No newline at end of file +TEST_CASE( "Test bivariate gaussian copula", "[prob_copula]" ) +{ + auto dist1 = Seldon::power_law_distribution( 0.02, 2.5 ); + auto dist2 = Seldon::truncated_normal_distribution( 1.0, 0.75, 0.2 ); + auto copula = Seldon::bivariate_gaussian_copula( 0.5, dist1, dist2 ); + write_results_to_file( 10000, copula, "gaussian_copula.txt" ); +}