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 <[email protected]>
seldon-code · Mar 25, 2024 · 376bb90 · 376bb90
1 parent ea7c900
commit 376bb90
Show file tree

Hide file tree

Showing 2 changed files with 14 additions and 32 deletions.
diff --git a/include/util/math.hpp b/include/util/math.hpp
@@ -243,11 +243,11 @@ class bivariate_gaussian_copula
 {
 private:
     ScalarT covariance;
-    bivariate_normal_distribution<ScalarT> bivariate_normal_dist{};
+    bivariate_normal_distribution<ScalarT> biv_normal_dist{};
     // std::normal_distribution<ScalarT> 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,22 +258,21 @@ 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 )
     {
     }
 
     template<typename Generator>
     std::array<ScalarT, 2> 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<ScalarT, 2> z_unit = { cum_gauss( z[0] ), cum_gauss( z[1] ) };
+        std::array<ScalarT, 2> z_unit = { cdf_gauss( z[0] ), cdf_gauss( z[1] ) };
         // 3. Apply inverse transform sampling
-        std::array<ScalarT, 2> res
-            = { dist1.inverse_cumulative_probability( z_unit[0] ), dist2.inverse_cumulative_probability( z_unit[1] ) };
+        std::array<ScalarT, 2> 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
@@ -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<ActivityDrivenModel::AgentT>( network_file );
-
-//     std::vector<double> opinions_expected    = { 2.1127107987061544, 0.8088982488089491, -0.8802809369462433 };
-//     std::vector<double> activities_expected  = { 0.044554683389757696, 0.015813166022685163, 0.015863953902810535 };
-//     std::vector<double> 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 ) );
-//     }
-// }
+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" );
+}