From 77878730f85c1c8beacf6b6dbc0f04822657d31e Mon Sep 17 00:00:00 2001 From: James Yang Date: Tue, 14 Apr 2020 18:40:18 -0400 Subject: [PATCH 01/35] Reorganize include dir --- CMakeLists.txt | 18 ++++++++++++++++-- autoppl/CMakeLists.txt | 15 --------------- include/autoppl/autoppl | 2 ++ .../autoppl/dummy/fib.hpp | 2 +- test/test1.cpp | 2 +- 5 files changed, 20 insertions(+), 19 deletions(-) delete mode 100644 autoppl/CMakeLists.txt create mode 100644 include/autoppl/autoppl rename autoppl/include/autoppl.hpp => include/autoppl/dummy/fib.hpp (84%) diff --git a/CMakeLists.txt b/CMakeLists.txt index 1d10b546..06f4fa53 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -10,7 +10,6 @@ option(AUTOPPL_ENABLE_BENCHMARK "Enable benchmarks to be built." OFF) option(AUTOPPL_ENABLE_TEST_COVERAGE "Build with test coverage (AUTOPPL_ENABLE_TEST must be ON)" OFF) # Automate the choosing of config -# if CMAKE_BUILD_TYPE not defined if (NOT CMAKE_BUILD_TYPE) # if binary directory ends with "release", use release mode if (${PROJECT_BINARY_DIR} MATCHES "release$") @@ -22,6 +21,16 @@ if (NOT CMAKE_BUILD_TYPE) endif() message(STATUS "Compiling in ${CMAKE_BUILD_TYPE} mode") +# Add this library as interface (header-only) +add_library(${PROJECT_NAME} INTERFACE) + +target_include_directories(${PROJECT_NAME} + INTERFACE $ + $) + +# Set C++17 standard for project target +target_compile_features(${PROJECT_NAME} INTERFACE cxx_std_17) + # Configure tests if (AUTOPPL_ENABLE_TEST) include(CTest) # enable memcheck @@ -35,6 +44,11 @@ if (AUTOPPL_ENABLE_TEST) add_subdirectory(${PROJECT_SOURCE_DIR}/test ${PROJECT_BINARY_DIR}/test) endif() +# TODO: add src dir if needed +#set(AUTOPPL_SOURCE_DIR ${CMAKE_CURRENT_LIST_DIR}/src) +#file(GLOB_RECURSE AUTOPPL_SOURCE_FILES RELATIVE src LIST_DIRECTORIES false *.cpp) +#set(AUTOPPL_SOURCE_FILES ${AUTOPPL_SOURCE_DIR}/autoppl.cpp) +#set(AUTOPPL_HEADER_FILES ${AUTOPPL_INCLUDE_DIR}/autoppl.h) + # Add subdirectories -add_subdirectory(${PROJECT_SOURCE_DIR}/autoppl ${PROJECT_BINARY_DIR}/autoppl) add_subdirectory(${PROJECT_SOURCE_DIR}/lib ${PROJECT_BINARY_DIR}/lib) diff --git a/autoppl/CMakeLists.txt b/autoppl/CMakeLists.txt deleted file mode 100644 index b50a56a2..00000000 --- a/autoppl/CMakeLists.txt +++ /dev/null @@ -1,15 +0,0 @@ -# Add this library as interface (header-only) -add_library(${PROJECT_NAME} INTERFACE) - -target_include_directories(${PROJECT_NAME} - INTERFACE $ - $) - -# Set C++17 standard for project target -target_compile_features(${PROJECT_NAME} INTERFACE cxx_std_17) - -#set(AUTOPPL_SOURCE_DIR ${CMAKE_CURRENT_LIST_DIR}/src) -#file(GLOB_RECURSE AUTOPPL_SOURCE_FILES RELATIVE src LIST_DIRECTORIES false *.cpp) - -#set(AUTOPPL_SOURCE_FILES ${AUTOPPL_SOURCE_DIR}/autoppl.cpp) -#set(AUTOPPL_HEADER_FILES ${AUTOPPL_INCLUDE_DIR}/autoppl.h) diff --git a/include/autoppl/autoppl b/include/autoppl/autoppl new file mode 100644 index 00000000..311a648c --- /dev/null +++ b/include/autoppl/autoppl @@ -0,0 +1,2 @@ +#pragma once +#include "bits/fib.hpp" diff --git a/autoppl/include/autoppl.hpp b/include/autoppl/dummy/fib.hpp similarity index 84% rename from autoppl/include/autoppl.hpp rename to include/autoppl/dummy/fib.hpp index 729eb861..1625a9bd 100644 --- a/autoppl/include/autoppl.hpp +++ b/include/autoppl/dummy/fib.hpp @@ -8,4 +8,4 @@ inline int fib(int n) else return fib(n-1) + fib(n-2); } -} // namespace autoppl +} // namespace ppl diff --git a/test/test1.cpp b/test/test1.cpp index 876c350c..d24ea926 100644 --- a/test/test1.cpp +++ b/test/test1.cpp @@ -1,5 +1,5 @@ #include "gtest/gtest.h" -#include "autoppl.hpp" +#include namespace { From a7a437bdcbeadaac123a87f8c2767f003d0c4290 Mon Sep 17 00:00:00 2001 From: James Yang Date: Tue, 14 Apr 2020 21:41:09 -0400 Subject: [PATCH 02/35] Finish proof of concept --- include/autoppl/autoppl | 2 +- include/autoppl/dummy/fib.hpp | 11 ---- include/autoppl/expression/core.hpp | 17 ++++++ include/autoppl/expression/dist_utils.hpp | 5 ++ include/autoppl/expression/model.hpp | 68 +++++++++++++++++++++++ include/autoppl/expression/rv_tag.hpp | 38 +++++++++++++ include/autoppl/expression/traits.hpp | 29 ++++++++++ include/autoppl/expression/uniform.hpp | 44 +++++++++++++++ test/CMakeLists.txt | 16 +++--- test/expression/core_unittest.cpp | 50 +++++++++++++++++ test/test1.cpp | 11 ---- 11 files changed, 260 insertions(+), 31 deletions(-) delete mode 100644 include/autoppl/dummy/fib.hpp create mode 100644 include/autoppl/expression/core.hpp create mode 100644 include/autoppl/expression/dist_utils.hpp create mode 100644 include/autoppl/expression/model.hpp create mode 100644 include/autoppl/expression/rv_tag.hpp create mode 100644 include/autoppl/expression/traits.hpp create mode 100644 include/autoppl/expression/uniform.hpp create mode 100644 test/expression/core_unittest.cpp delete mode 100644 test/test1.cpp diff --git a/include/autoppl/autoppl b/include/autoppl/autoppl index 311a648c..1532ee7e 100644 --- a/include/autoppl/autoppl +++ b/include/autoppl/autoppl @@ -1,2 +1,2 @@ #pragma once -#include "bits/fib.hpp" +// TODO: export all headers later! diff --git a/include/autoppl/dummy/fib.hpp b/include/autoppl/dummy/fib.hpp deleted file mode 100644 index 1625a9bd..00000000 --- a/include/autoppl/dummy/fib.hpp +++ /dev/null @@ -1,11 +0,0 @@ -#pragma once - -namespace ppl { - -inline int fib(int n) -{ - if (n <= 1) return 1; - else return fib(n-1) + fib(n-2); -} - -} // namespace ppl diff --git a/include/autoppl/expression/core.hpp b/include/autoppl/expression/core.hpp new file mode 100644 index 00000000..f44dde3e --- /dev/null +++ b/include/autoppl/expression/core.hpp @@ -0,0 +1,17 @@ +#pragma once +#include +#include +#include + +namespace ppl { + +// Builds an EqNode to associate tag with dist. +// Ex. x |= uniform(0,1) +template +constexpr inline auto operator|=(rv_tag& tag, + DistType dist) +{ + return EqNode, DistType>(tag, dist); +} + +} // namespace ppl diff --git a/include/autoppl/expression/dist_utils.hpp b/include/autoppl/expression/dist_utils.hpp new file mode 100644 index 00000000..df05f7ce --- /dev/null +++ b/include/autoppl/expression/dist_utils.hpp @@ -0,0 +1,5 @@ +#pragma once + +namespace ppl { + +} // namespace ppl diff --git a/include/autoppl/expression/model.hpp b/include/autoppl/expression/model.hpp new file mode 100644 index 00000000..265b9455 --- /dev/null +++ b/include/autoppl/expression/model.hpp @@ -0,0 +1,68 @@ +#pragma once +#include +#include + +namespace ppl { +namespace details { + +/* + * The possible states for a node expression. + * By default, all nodes should be a parameter. + */ +enum class node_state : bool { + data, + parameter +}; + +} // namespace details + +/* + * This class represents a "node" in the model expression + * that relates a tag with a distribution. + */ +template +struct EqNode +{ + using tag_t = TagType; + using tag_cref_t = std::reference_wrapper; + using value_t = typename tag_traits::value_t; + using pointer_t = typename tag_traits::pointer_t; + + using state_t = details::node_state; + + using dist_t = DistType; + using gen_value_t = typename tag_traits::value_t; + using dist_value_t = typename dist_traits::dist_value_t; + + EqNode(const tag_t& tag, dist_t dist) noexcept + : value_{} + , tag_{} + , state_{state_t::parameter} + , tag_cref_{tag} + , dist_{dist} + {} + + /* + * Compute pdf of dist_ with value x. + */ + dist_value_t pdf(gen_value_t x) const + { return dist_.pdf(x); } + + /* + * Compute log-pdf of dist_ with value x. + */ + dist_value_t log_pdf(gen_value_t x) const + { return dist_.log_pdf(x); } + +private: + // cache optimization + value_t value_; // value to store during computation + tag_t tag_; // tag to manage the storage: get/put values + + state_t state_; // state to determine if data or param + tag_cref_t tag_cref_; // (const) reference of the tag since + // storage of tag may be bound right before some computation + dist_t dist_; // distribution associated with tag +}; + +} // namespace ppl diff --git a/include/autoppl/expression/rv_tag.hpp b/include/autoppl/expression/rv_tag.hpp new file mode 100644 index 00000000..1026c957 --- /dev/null +++ b/include/autoppl/expression/rv_tag.hpp @@ -0,0 +1,38 @@ +#pragma once + +namespace ppl { + +/* + * rv_tag is a light-weight structure that represents a univariate random variable. + * It acts as an intermediate layer of communication between + * a model expression and the users, who must supply storage of values associated with this tag. + */ +template +struct rv_tag +{ + using value_t = ValueType; + using pointer_t = value_t*; + + rv_tag(pointer_t storage_ptr = nullptr) noexcept + : storage_ptr_{storage_ptr} + {} + + /* + * Binds storage pointer to storage_ptr. + * Resets pointer to original source to be storage_ptr. + */ + void bind_storage(pointer_t storage_ptr) + { + storage_ptr_ = storage_ptr; + } + +private: + pointer_t storage_ptr_; // points to beginning of storage + // storage is assumed to be contiguous +}; + +// Useful aliases +using cont_rv_tag = rv_tag; // continuous RV tag +using disc_rv_tag = rv_tag; // discrete RV tag + +} // namespace ppl diff --git a/include/autoppl/expression/traits.hpp b/include/autoppl/expression/traits.hpp new file mode 100644 index 00000000..8680ea83 --- /dev/null +++ b/include/autoppl/expression/traits.hpp @@ -0,0 +1,29 @@ +#pragma once + +namespace ppl { + +/* + * This class lists all member aliases that + * any tag types should have, in particular, rv_tag. + * Users should rely on this class to grab member aliases. + */ +template +struct tag_traits +{ + using value_t = typename TagType::value_t; + using pointer_t = typename TagType::pointer_t; +}; + +/* + * This class lists all member aliases that + * any distribution types should have. + * Users should rely on this class to grab member aliases. + */ +template +struct dist_traits +{ + using value_t = typename DistType::value_t; + using dist_value_t = typename DistType::dist_value_t; +}; + +} // namespace ppl diff --git a/include/autoppl/expression/uniform.hpp b/include/autoppl/expression/uniform.hpp new file mode 100644 index 00000000..3d1d6f50 --- /dev/null +++ b/include/autoppl/expression/uniform.hpp @@ -0,0 +1,44 @@ +#pragma once +#include +#include +#include +#include + +namespace ppl { + +struct uniform +{ + using value_t = double; + using dist_value_t = double; + + uniform(value_t min, value_t max) + : min_{min}, max_{max} + { + assert(min_ < max_); + } + + // TODO: tag this class as "TriviallySamplable" + template + value_t sample(GeneratorType& gen) const + { + std::uniform_real_distribution dist(min_, max_); + return dist(gen); + } + + dist_value_t pdf(value_t x) const + { + return (min_ < x && x < max_) ? 1./(max_ - min_) : 0; + } + + dist_value_t log_pdf(value_t x) const + { + return (min_ < x && x < max_) ? + -std::log(max_ - min_) : + std::numeric_limits::lowest(); + } + +private: + value_t min_, max_; +}; + +} // namespace ppl diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 7a21c768..8b1247e0 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -6,16 +6,16 @@ if (AUTOPPL_ENABLE_TEST_COVERAGE) endif() ###################################################### -# Dummy Test +# Expression Test ###################################################### -add_executable(dummy_unittest - test1.cpp +add_executable(expression_unittest + ${CMAKE_CURRENT_SOURCE_DIR}/expression/core_unittest.cpp ) -target_compile_options(dummy_unittest PRIVATE -g -Wall -Werror -Wextra) -target_include_directories(dummy_unittest PRIVATE ${GTEST_DIR}/include) +target_compile_options(expression_unittest PRIVATE -g -Wall -Werror -Wextra) +target_include_directories(expression_unittest PRIVATE ${GTEST_DIR}/include) if (AUTOPPL_ENABLE_TEST_COVERAGE) - target_link_libraries(dummy_unittest gcov) + target_link_libraries(expression_unittest gcov) endif() -target_link_libraries(dummy_unittest gtest_main pthread ${PROJECT_NAME}) -add_test(dummy_unittest dummy_unittest) +target_link_libraries(expression_unittest gtest_main pthread ${PROJECT_NAME}) +add_test(expression_unittest expression_unittest) diff --git a/test/expression/core_unittest.cpp b/test/expression/core_unittest.cpp new file mode 100644 index 00000000..5d7789a5 --- /dev/null +++ b/test/expression/core_unittest.cpp @@ -0,0 +1,50 @@ +#include "gtest/gtest.h" +#include + +namespace ppl { + +/* + * Fixture for testing rv_tag with distribution. + */ +struct tag_dist_fixture : ::testing::Test +{ +protected: + rv_tag x; +}; + +TEST_F(tag_dist_fixture, pdf_valid) +{ + auto model = (x |= uniform(0,1)); + EXPECT_EQ(model.pdf(0.000001), 1); + EXPECT_EQ(model.pdf(0.5), 1); + EXPECT_EQ(model.pdf(0.999999), 1); +} + +TEST_F(tag_dist_fixture, pdf_invalid) +{ + auto model = (x |= uniform(0,1)); + EXPECT_EQ(model.pdf(-69), 0); + EXPECT_EQ(model.pdf(-0.00001), 0); + EXPECT_EQ(model.pdf(1.00001), 0); + EXPECT_EQ(model.pdf(69), 0); +} + +TEST_F(tag_dist_fixture, log_pdf_valid) +{ + auto model = (x |= uniform(0,1)); + EXPECT_EQ(model.log_pdf(0.000001), 0); + EXPECT_EQ(model.log_pdf(0.5), 0); + EXPECT_EQ(model.log_pdf(0.999999), 0); +} + +TEST_F(tag_dist_fixture, log_pdf_invalid) +{ + auto model = (x |= uniform(0,1)); + double expected = std::numeric_limits::lowest(); + EXPECT_EQ(model.log_pdf(-69), expected); + EXPECT_EQ(model.log_pdf(-0.00001), expected); + EXPECT_EQ(model.log_pdf(1.00001), expected); + EXPECT_EQ(model.log_pdf(69), expected); +} + +} // namespace ppl diff --git a/test/test1.cpp b/test/test1.cpp deleted file mode 100644 index d24ea926..00000000 --- a/test/test1.cpp +++ /dev/null @@ -1,11 +0,0 @@ -#include "gtest/gtest.h" -#include - -namespace { - -TEST(blaTest, test1) { - int n = ppl::fib(10); - EXPECT_EQ(n, 89); -} - -} From e041f6b118e3b92f28709cfe820c8b1e9a49a048 Mon Sep 17 00:00:00 2001 From: James Yang Date: Tue, 14 Apr 2020 21:42:57 -0400 Subject: [PATCH 03/35] Remove unnecessary file --- include/autoppl/expression/dist_utils.hpp | 5 ----- 1 file changed, 5 deletions(-) delete mode 100644 include/autoppl/expression/dist_utils.hpp diff --git a/include/autoppl/expression/dist_utils.hpp b/include/autoppl/expression/dist_utils.hpp deleted file mode 100644 index df05f7ce..00000000 --- a/include/autoppl/expression/dist_utils.hpp +++ /dev/null @@ -1,5 +0,0 @@ -#pragma once - -namespace ppl { - -} // namespace ppl From 09fe7df8bf9b8f0abe5f5694428065027802b624 Mon Sep 17 00:00:00 2001 From: James Yang Date: Tue, 14 Apr 2020 21:47:26 -0400 Subject: [PATCH 04/35] Remove outdated comment about bind_storage --- include/autoppl/expression/rv_tag.hpp | 5 +---- 1 file changed, 1 insertion(+), 4 deletions(-) diff --git a/include/autoppl/expression/rv_tag.hpp b/include/autoppl/expression/rv_tag.hpp index 1026c957..b245656d 100644 --- a/include/autoppl/expression/rv_tag.hpp +++ b/include/autoppl/expression/rv_tag.hpp @@ -19,12 +19,9 @@ struct rv_tag /* * Binds storage pointer to storage_ptr. - * Resets pointer to original source to be storage_ptr. */ void bind_storage(pointer_t storage_ptr) - { - storage_ptr_ = storage_ptr; - } + { storage_ptr_ = storage_ptr; } private: pointer_t storage_ptr_; // points to beginning of storage From d4692b32b1991194095c647dedfdd5e5adb87114 Mon Sep 17 00:00:00 2001 From: James Yang Date: Tue, 14 Apr 2020 21:52:54 -0400 Subject: [PATCH 05/35] Update travis coveralls include dir --- .travis.yml | 2 +- test/expression/core_unittest.cpp | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/.travis.yml b/.travis.yml index e7130246..54bc91a1 100644 --- a/.travis.yml +++ b/.travis.yml @@ -38,7 +38,7 @@ jobs: coveralls --root ../../ --build-root ./ - --include autoppl/include + --include include --exclude lib --gcov 'gcov-7' --gcov-options '\-lp' diff --git a/test/expression/core_unittest.cpp b/test/expression/core_unittest.cpp index 5d7789a5..f39debdc 100644 --- a/test/expression/core_unittest.cpp +++ b/test/expression/core_unittest.cpp @@ -4,7 +4,7 @@ namespace ppl { /* - * Fixture for testing rv_tag with distribution. + * Fixture for testing one rv_tag with distribution. */ struct tag_dist_fixture : ::testing::Test { From ebe01dc2a94c3501a5a81f78a7078c49c08cb052 Mon Sep 17 00:00:00 2001 From: James Yang Date: Wed, 15 Apr 2020 00:02:45 -0400 Subject: [PATCH 06/35] Added new features to allow for glueing expressions --- include/autoppl/expression/core.hpp | 24 +++++- include/autoppl/expression/model.hpp | 106 +++++++++++++++++-------- include/autoppl/expression/rv_tag.hpp | 53 +++++++++++-- include/autoppl/expression/traits.hpp | 19 +++-- include/autoppl/expression/uniform.hpp | 8 +- test/expression/core_unittest.cpp | 98 ++++++++++++++++++----- 6 files changed, 233 insertions(+), 75 deletions(-) diff --git a/include/autoppl/expression/core.hpp b/include/autoppl/expression/core.hpp index f44dde3e..4fcc0e1f 100644 --- a/include/autoppl/expression/core.hpp +++ b/include/autoppl/expression/core.hpp @@ -5,13 +5,29 @@ namespace ppl { -// Builds an EqNode to associate tag with dist. -// Ex. x |= uniform(0,1) +// TODO: all these template parameters should be constrained +// with concepts! + +/* + * Builds an EqNode to associate tag with dist. + * Ex. x |= uniform(0,1) + */ template -constexpr inline auto operator|=(rv_tag& tag, - DistType dist) +constexpr inline auto operator|=(const rv_tag& tag, + const DistType& dist) { return EqNode, DistType>(tag, dist); } +/* + * Builds a GlueNode to "glue" the left expression with the right. + * Ex. (x |= uniform(0,1), y |= uniform(0, 2)) + */ +template +constexpr inline auto operator,(const LHSNodeType& lhs, + const RHSNodeType& rhs) +{ + return GlueNode(lhs, rhs); +} + } // namespace ppl diff --git a/include/autoppl/expression/model.hpp b/include/autoppl/expression/model.hpp index 265b9455..5a4714d7 100644 --- a/include/autoppl/expression/model.hpp +++ b/include/autoppl/expression/model.hpp @@ -1,20 +1,9 @@ #pragma once +#include #include #include namespace ppl { -namespace details { - -/* - * The possible states for a node expression. - * By default, all nodes should be a parameter. - */ -enum class node_state : bool { - data, - parameter -}; - -} // namespace details /* * This class represents a "node" in the model expression @@ -24,45 +13,92 @@ template struct EqNode { using tag_t = TagType; - using tag_cref_t = std::reference_wrapper; - using value_t = typename tag_traits::value_t; - using pointer_t = typename tag_traits::pointer_t; - - using state_t = details::node_state; - using dist_t = DistType; - using gen_value_t = typename tag_traits::value_t; using dist_value_t = typename dist_traits::dist_value_t; - EqNode(const tag_t& tag, dist_t dist) noexcept - : value_{} - , tag_{} - , state_{state_t::parameter} + EqNode(const tag_t& tag, + const dist_t& dist) noexcept + : tag_{} , tag_cref_{tag} , dist_{dist} {} /* - * Compute pdf of dist_ with value x. + * Updates underlying tag by copying from referenced tag. + * This must be called before the model is used, if + * any members of the referenced tag changed. */ - dist_value_t pdf(gen_value_t x) const - { return dist_.pdf(x); } + void update() + { tag_ = tag_cref_; } + + /* + * Compute pdf of dist_ with value_. + * Assumes that value_ has been assigned with a proper value. + */ + dist_value_t pdf() const + { return dist_.pdf(tag_.get_value()); } /* * Compute log-pdf of dist_ with value x. + * Assumes that value_ has been assigned with a proper value. */ - dist_value_t log_pdf(gen_value_t x) const - { return dist_.log_pdf(x); } + dist_value_t log_pdf() const + { return dist_.log_pdf(tag_.get_value()); } private: - // cache optimization - value_t value_; // value to store during computation - tag_t tag_; // tag to manage the storage: get/put values - - state_t state_; // state to determine if data or param - tag_cref_t tag_cref_; // (const) reference of the tag since - // storage of tag may be bound right before some computation + using tag_cref_t = std::reference_wrapper; + + tag_t tag_; // cache optimization + tag_cref_t tag_cref_; // (const) reference of the tag since any configuration + // may be changed until right before update dist_t dist_; // distribution associated with tag }; +/* + * This class represents a "node" in a model expression that + * "glues" two sub-model expressions. + */ +template +struct GlueNode +{ + using left_node_t = LHSNodeType; + using right_node_t = RHSNodeType; + using dist_value_t = std::common_type_t< + typename node_traits::dist_value_t, + typename node_traits::dist_value_t + >; + + GlueNode(const left_node_t& lhs, + const right_node_t& rhs) + : left_node_{lhs} + , right_node_{rhs} + {} + + /* + * Updates left node first then right node by calling "update" on each. + * This must be called before the model is used, + * if either left or right must be updated. + */ + void update() + { left_node_.update(); right_node_.update(); } + + /* + * Computes left node joint pdf then right node joint pdf + * and returns the product of the two. + */ + dist_value_t pdf() const + { return left_node_.pdf() * right_node_.pdf(); } + + /* + * Computes left node joint log-pdf then right node joint log-pdf + * and returns the sum of the two. + */ + dist_value_t log_pdf() const + { return left_node_.log_pdf() + right_node_.log_pdf(); } + +private: + left_node_t left_node_; + right_node_t right_node_; +}; + } // namespace ppl diff --git a/include/autoppl/expression/rv_tag.hpp b/include/autoppl/expression/rv_tag.hpp index b245656d..00ccaa65 100644 --- a/include/autoppl/expression/rv_tag.hpp +++ b/include/autoppl/expression/rv_tag.hpp @@ -2,6 +2,16 @@ namespace ppl { +/* + * The possible states for a tag. + * By default, all tags should be considered as a parameter. + * TODO: maybe move in a different file? + */ +enum class tag_state : bool { + data, + parameter +}; + /* * rv_tag is a light-weight structure that represents a univariate random variable. * It acts as an intermediate layer of communication between @@ -12,20 +22,51 @@ struct rv_tag { using value_t = ValueType; using pointer_t = value_t*; + using const_pointer_t = const value_t*; + using state_t = tag_state; + + // constructors + rv_tag(value_t value, + pointer_t storage_ptr) noexcept + : value_{value} + , storage_ptr_{storage_ptr} + , state_{state_t::parameter} + {} - rv_tag(pointer_t storage_ptr = nullptr) noexcept - : storage_ptr_{storage_ptr} + rv_tag(pointer_t storage_ptr) noexcept + : rv_tag(0, storage_ptr) {} + rv_tag() noexcept + : rv_tag(0, nullptr) + {} + + void set_value(value_t value) { value_ = value; } + value_t get_value() const { return value_; } + + void set_storage(pointer_t storage_ptr) { storage_ptr_ = storage_ptr; } + pointer_t get_storage() { return storage_ptr_; } + const_pointer_t get_storage() const { return storage_ptr_; } + + void set_state(state_t state) { state_ = state; } + state_t get_state() const { return state_; } + /* - * Binds storage pointer to storage_ptr. + * Sets underlying value to "value". + * Additionally modifies the tag to be considered as data. + * Equivalent to calling set_value(value) then set_state(state). */ - void bind_storage(pointer_t storage_ptr) - { storage_ptr_ = storage_ptr; } + void observe(value_t value) + { + set_value(value); + set_state(state_t::data); + } private: - pointer_t storage_ptr_; // points to beginning of storage + value_t value_; // store value associated with tag + pointer_t storage_ptr_; // points to beginning of storage // storage is assumed to be contiguous + state_t state_; // state to determine if data or param }; // Useful aliases diff --git a/include/autoppl/expression/traits.hpp b/include/autoppl/expression/traits.hpp index 8680ea83..dce5efc6 100644 --- a/include/autoppl/expression/traits.hpp +++ b/include/autoppl/expression/traits.hpp @@ -3,22 +3,19 @@ namespace ppl { /* - * This class lists all member aliases that - * any tag types should have, in particular, rv_tag. - * Users should rely on this class to grab member aliases. + * The following classes list member aliases that + * any such template parameter types should have. + * Users should rely on these classes to grab member aliases. */ + template struct tag_traits { using value_t = typename TagType::value_t; using pointer_t = typename TagType::pointer_t; + using state_t = typename TagType::state_t; }; -/* - * This class lists all member aliases that - * any distribution types should have. - * Users should rely on this class to grab member aliases. - */ template struct dist_traits { @@ -26,4 +23,10 @@ struct dist_traits using dist_value_t = typename DistType::dist_value_t; }; +template +struct node_traits +{ + using dist_value_t = typename NodeType::dist_value_t; +}; + } // namespace ppl diff --git a/include/autoppl/expression/uniform.hpp b/include/autoppl/expression/uniform.hpp index 3d1d6f50..562d2c6b 100644 --- a/include/autoppl/expression/uniform.hpp +++ b/include/autoppl/expression/uniform.hpp @@ -6,6 +6,8 @@ namespace ppl { +// TODO: change name to UniformDist and make class template. +// uniform should be a function that creates this kind of object. struct uniform { using value_t = double; @@ -13,11 +15,9 @@ struct uniform uniform(value_t min, value_t max) : min_{min}, max_{max} - { - assert(min_ < max_); - } + { assert(min_ < max_); } - // TODO: tag this class as "TriviallySamplable" + // TODO: tag this class as "TriviallySamplable"? template value_t sample(GeneratorType& gen) const { diff --git a/test/expression/core_unittest.cpp b/test/expression/core_unittest.cpp index f39debdc..2ac93e7c 100644 --- a/test/expression/core_unittest.cpp +++ b/test/expression/core_unittest.cpp @@ -3,6 +3,10 @@ namespace ppl { +////////////////////////////////////////////////////// +// Model with one RV TESTS +////////////////////////////////////////////////////// + /* * Fixture for testing one rv_tag with distribution. */ @@ -10,41 +14,99 @@ struct tag_dist_fixture : ::testing::Test { protected: rv_tag x; + using model_t = std::decay_t; + model_t model = (x |= uniform(0,1)); + + void reconfigure(double val) + { + x.set_value(val); + model.update(); + } }; TEST_F(tag_dist_fixture, pdf_valid) { - auto model = (x |= uniform(0,1)); - EXPECT_EQ(model.pdf(0.000001), 1); - EXPECT_EQ(model.pdf(0.5), 1); - EXPECT_EQ(model.pdf(0.999999), 1); + reconfigure(0.000001); + EXPECT_EQ(model.pdf(), 1); + reconfigure(0.5); + EXPECT_EQ(model.pdf(), 1); + reconfigure(0.999999); + EXPECT_EQ(model.pdf(), 1); } TEST_F(tag_dist_fixture, pdf_invalid) { - auto model = (x |= uniform(0,1)); - EXPECT_EQ(model.pdf(-69), 0); - EXPECT_EQ(model.pdf(-0.00001), 0); - EXPECT_EQ(model.pdf(1.00001), 0); - EXPECT_EQ(model.pdf(69), 0); + reconfigure(-69); + EXPECT_EQ(model.pdf(), 0); + reconfigure(-0.000001); + EXPECT_EQ(model.pdf(), 0); + reconfigure(1.000001); + EXPECT_EQ(model.pdf(), 0); + reconfigure(69); + EXPECT_EQ(model.pdf(), 0); } TEST_F(tag_dist_fixture, log_pdf_valid) { - auto model = (x |= uniform(0,1)); - EXPECT_EQ(model.log_pdf(0.000001), 0); - EXPECT_EQ(model.log_pdf(0.5), 0); - EXPECT_EQ(model.log_pdf(0.999999), 0); + reconfigure(0.000001); + EXPECT_EQ(model.log_pdf(), 0); + reconfigure(0.5); + EXPECT_EQ(model.log_pdf(), 0); + reconfigure(0.999999); + EXPECT_EQ(model.log_pdf(), 0); } TEST_F(tag_dist_fixture, log_pdf_invalid) { - auto model = (x |= uniform(0,1)); double expected = std::numeric_limits::lowest(); - EXPECT_EQ(model.log_pdf(-69), expected); - EXPECT_EQ(model.log_pdf(-0.00001), expected); - EXPECT_EQ(model.log_pdf(1.00001), expected); - EXPECT_EQ(model.log_pdf(69), expected); + reconfigure(-69); + EXPECT_EQ(model.log_pdf(), expected); + reconfigure(-0.000001); + EXPECT_EQ(model.log_pdf(), expected); + reconfigure(1.000001); + EXPECT_EQ(model.log_pdf(), expected); + reconfigure(69); + EXPECT_EQ(model.log_pdf(), expected); +} + +////////////////////////////////////////////////////// +// Model with many RV (no dependencies) TESTS +////////////////////////////////////////////////////// + +struct many_tag_dist_fixture : ::testing::Test +{ +protected: + rv_tag x, y, z, w; +}; + +TEST_F(many_tag_dist_fixture, two_tags) +{ + auto model = ( + x |= uniform(0, 1), + y |= uniform(0, 2) + ); + + auto reconfigure = [&](double xv, double yv) { + x.set_value(xv); + y.set_value(yv); + model.update(); + }; + + // both within range + reconfigure(0.2, 1.8); + EXPECT_EQ(model.pdf(), 0.5); + + // first out of range + reconfigure(-0.0005, 0.99999); + EXPECT_EQ(model.pdf(), 0); + + // second out of range + reconfigure(0.5142, 2.0000123); + EXPECT_EQ(model.pdf(), 0); + + // both out of range + reconfigure(-10., 69.); + EXPECT_EQ(model.pdf(), 0); } } // namespace ppl From 98d88ca972ca69271e5e00fa234e01492af2e3b5 Mon Sep 17 00:00:00 2001 From: James Yang Date: Wed, 15 Apr 2020 00:34:20 -0400 Subject: [PATCH 07/35] Correctly mocked out distribution and tag types in testing --- include/autoppl/expression/core.hpp | 33 ----- include/autoppl/expression/model.hpp | 29 +++++ test/CMakeLists.txt | 2 +- test/expression/core_unittest.cpp | 112 ---------------- test/expression/model_unittest.cpp | 184 +++++++++++++++++++++++++++ 5 files changed, 214 insertions(+), 146 deletions(-) delete mode 100644 include/autoppl/expression/core.hpp delete mode 100644 test/expression/core_unittest.cpp create mode 100644 test/expression/model_unittest.cpp diff --git a/include/autoppl/expression/core.hpp b/include/autoppl/expression/core.hpp deleted file mode 100644 index 4fcc0e1f..00000000 --- a/include/autoppl/expression/core.hpp +++ /dev/null @@ -1,33 +0,0 @@ -#pragma once -#include -#include -#include - -namespace ppl { - -// TODO: all these template parameters should be constrained -// with concepts! - -/* - * Builds an EqNode to associate tag with dist. - * Ex. x |= uniform(0,1) - */ -template -constexpr inline auto operator|=(const rv_tag& tag, - const DistType& dist) -{ - return EqNode, DistType>(tag, dist); -} - -/* - * Builds a GlueNode to "glue" the left expression with the right. - * Ex. (x |= uniform(0,1), y |= uniform(0, 2)) - */ -template -constexpr inline auto operator,(const LHSNodeType& lhs, - const RHSNodeType& rhs) -{ - return GlueNode(lhs, rhs); -} - -} // namespace ppl diff --git a/include/autoppl/expression/model.hpp b/include/autoppl/expression/model.hpp index 5a4714d7..dabb7fc2 100644 --- a/include/autoppl/expression/model.hpp +++ b/include/autoppl/expression/model.hpp @@ -101,4 +101,33 @@ struct GlueNode right_node_t right_node_; }; +///////////////////////////////////////////////////////// +// Operator overloads +///////////////////////////////////////////////////////// + +// TODO: all these template parameters should be constrained +// with concepts! + +/* + * Builds an EqNode to associate tag with dist. + * Ex. x |= uniform(0,1) + */ +template +constexpr inline auto operator|=(const TagType& tag, + const DistType& dist) +{ + return EqNode(tag, dist); +} + +/* + * Builds a GlueNode to "glue" the left expression with the right. + * Ex. (x |= uniform(0,1), y |= uniform(0, 2)) + */ +template +constexpr inline auto operator,(const LHSNodeType& lhs, + const RHSNodeType& rhs) +{ + return GlueNode(lhs, rhs); +} + } // namespace ppl diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 8b1247e0..ef3cd0e2 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -10,7 +10,7 @@ endif() ###################################################### add_executable(expression_unittest - ${CMAKE_CURRENT_SOURCE_DIR}/expression/core_unittest.cpp + ${CMAKE_CURRENT_SOURCE_DIR}/expression/model_unittest.cpp ) target_compile_options(expression_unittest PRIVATE -g -Wall -Werror -Wextra) target_include_directories(expression_unittest PRIVATE ${GTEST_DIR}/include) diff --git a/test/expression/core_unittest.cpp b/test/expression/core_unittest.cpp deleted file mode 100644 index 2ac93e7c..00000000 --- a/test/expression/core_unittest.cpp +++ /dev/null @@ -1,112 +0,0 @@ -#include "gtest/gtest.h" -#include - -namespace ppl { - -////////////////////////////////////////////////////// -// Model with one RV TESTS -////////////////////////////////////////////////////// - -/* - * Fixture for testing one rv_tag with distribution. - */ -struct tag_dist_fixture : ::testing::Test -{ -protected: - rv_tag x; - using model_t = std::decay_t; - model_t model = (x |= uniform(0,1)); - - void reconfigure(double val) - { - x.set_value(val); - model.update(); - } -}; - -TEST_F(tag_dist_fixture, pdf_valid) -{ - reconfigure(0.000001); - EXPECT_EQ(model.pdf(), 1); - reconfigure(0.5); - EXPECT_EQ(model.pdf(), 1); - reconfigure(0.999999); - EXPECT_EQ(model.pdf(), 1); -} - -TEST_F(tag_dist_fixture, pdf_invalid) -{ - reconfigure(-69); - EXPECT_EQ(model.pdf(), 0); - reconfigure(-0.000001); - EXPECT_EQ(model.pdf(), 0); - reconfigure(1.000001); - EXPECT_EQ(model.pdf(), 0); - reconfigure(69); - EXPECT_EQ(model.pdf(), 0); -} - -TEST_F(tag_dist_fixture, log_pdf_valid) -{ - reconfigure(0.000001); - EXPECT_EQ(model.log_pdf(), 0); - reconfigure(0.5); - EXPECT_EQ(model.log_pdf(), 0); - reconfigure(0.999999); - EXPECT_EQ(model.log_pdf(), 0); -} - -TEST_F(tag_dist_fixture, log_pdf_invalid) -{ - double expected = std::numeric_limits::lowest(); - reconfigure(-69); - EXPECT_EQ(model.log_pdf(), expected); - reconfigure(-0.000001); - EXPECT_EQ(model.log_pdf(), expected); - reconfigure(1.000001); - EXPECT_EQ(model.log_pdf(), expected); - reconfigure(69); - EXPECT_EQ(model.log_pdf(), expected); -} - -////////////////////////////////////////////////////// -// Model with many RV (no dependencies) TESTS -////////////////////////////////////////////////////// - -struct many_tag_dist_fixture : ::testing::Test -{ -protected: - rv_tag x, y, z, w; -}; - -TEST_F(many_tag_dist_fixture, two_tags) -{ - auto model = ( - x |= uniform(0, 1), - y |= uniform(0, 2) - ); - - auto reconfigure = [&](double xv, double yv) { - x.set_value(xv); - y.set_value(yv); - model.update(); - }; - - // both within range - reconfigure(0.2, 1.8); - EXPECT_EQ(model.pdf(), 0.5); - - // first out of range - reconfigure(-0.0005, 0.99999); - EXPECT_EQ(model.pdf(), 0); - - // second out of range - reconfigure(0.5142, 2.0000123); - EXPECT_EQ(model.pdf(), 0); - - // both out of range - reconfigure(-10., 69.); - EXPECT_EQ(model.pdf(), 0); -} - -} // namespace ppl diff --git a/test/expression/model_unittest.cpp b/test/expression/model_unittest.cpp new file mode 100644 index 00000000..62b2c05c --- /dev/null +++ b/test/expression/model_unittest.cpp @@ -0,0 +1,184 @@ +#include +#include "gtest/gtest.h" +#include + +namespace ppl { + +////////////////////////////////////////////////////// +// Model with one RV TESTS +////////////////////////////////////////////////////// + +/* + * Mock tag object for testing purposes. + * Must meet some of the requirements of actual tag types. + */ +struct MockTag +{ + using value_t = double; + using pointer_t = double*; + using state_t = void; + + void set_value(double val) { value_ = val; } + double get_value() const { return value_; } + +private: + double value_; +}; + +/* + * Mock distribution object for testing purposes. + * Must meet some of the requirements of actual distribution types. + */ +struct MockDist +{ + using value_t = double; + using dist_value_t = double; + + // Mock pdf - identity function. + double pdf(double x) const + { return x; } + + // Mock log_pdf - log(pdf(x)). + double log_pdf(double x) const + { return std::log(pdf(x)); } +}; + +/* + * Fixture for testing one rv_tag with distribution. + */ +struct tag_dist_fixture : ::testing::Test +{ +protected: + MockTag x; + using model_t = std::decay_t; + model_t model = (x |= MockDist()); + double val; + + void reconfigure(double val) + { + x.set_value(val); + model.update(); + } +}; + +TEST_F(tag_dist_fixture, pdf_valid) +{ + // MockDist pdf is identity function + // so we may simply compare model.pdf() with val. + + val = 0.000001; + reconfigure(val); + EXPECT_EQ(model.pdf(), val); + + val = 0.5; + reconfigure(val); + EXPECT_EQ(model.pdf(), val); + + val = 0.999999; + reconfigure(val); + EXPECT_EQ(model.pdf(), val); +} + +TEST_F(tag_dist_fixture, pdf_invalid) +{ + val = 0.000004123; + reconfigure(val); + EXPECT_EQ(model.pdf(), val); + + val = 0.55555; + reconfigure(val); + EXPECT_EQ(model.pdf(), val); + + val = 5.234424231; + reconfigure(val); + EXPECT_EQ(model.pdf(), val); + + val = 69; + reconfigure(val); + EXPECT_EQ(model.pdf(), val); +} + +TEST_F(tag_dist_fixture, log_pdf_valid) +{ + val = 0.000001; + reconfigure(val); + EXPECT_EQ(model.log_pdf(), std::log(val)); + + val = 0.5; + reconfigure(val); + EXPECT_EQ(model.log_pdf(), std::log(val)); + + val = 0.999999; + reconfigure(val); + EXPECT_EQ(model.log_pdf(), std::log(val)); +} + +TEST_F(tag_dist_fixture, log_pdf_invalid) +{ + val = 0.000004123; + reconfigure(val); + EXPECT_EQ(model.log_pdf(), std::log(val)); + + val = 0.55555; + reconfigure(val); + EXPECT_EQ(model.log_pdf(), std::log(val)); + + val = 5.234424231; + reconfigure(val); + EXPECT_EQ(model.log_pdf(), std::log(val)); + + val = 69; + reconfigure(val); + EXPECT_EQ(model.log_pdf(), std::log(val)); +} + +////////////////////////////////////////////////////// +// Model with many RV (no dependencies) TESTS +////////////////////////////////////////////////////// + +struct many_tag_dist_fixture : ::testing::Test +{ +protected: + MockTag x, y, z, w; + double xv, yv, zv, wv; +}; + +TEST_F(many_tag_dist_fixture, two_tags) +{ + auto model = ( + x |= MockDist(), + y |= MockDist() + ); + + auto reconfigure = [&](double xv, double yv) { + x.set_value(xv); + y.set_value(yv); + model.update(); + }; + + // both within range + xv = 0.2; yv = 1.8; + reconfigure(xv, yv); + EXPECT_EQ(model.pdf(), xv * yv); + EXPECT_EQ(model.log_pdf(), std::log(xv) + std::log(yv)); + + // first out of range + xv = 0.000542431; yv = 0.99999; + reconfigure(xv, yv); + EXPECT_EQ(model.pdf(), xv * yv); + EXPECT_EQ(model.log_pdf(), std::log(xv) + std::log(yv)); + + // second out of range + xv = 0.5142; yv = 2.0000123; + reconfigure(xv, yv); + EXPECT_EQ(model.pdf(), xv * yv); + EXPECT_EQ(model.log_pdf(), std::log(xv) + std::log(yv)); + + // both out of range + xv = 14.3; yv = 69.; + reconfigure(xv, yv); + EXPECT_EQ(model.pdf(), xv * yv); + EXPECT_EQ(model.log_pdf(), std::log(xv) + std::log(yv)); +} + +} // namespace ppl From c77e23de995cbfc652fd64cfd1c4d974b5cfbfa8 Mon Sep 17 00:00:00 2001 From: James Yang Date: Wed, 15 Apr 2020 01:11:51 -0400 Subject: [PATCH 08/35] Added one more test --- include/autoppl/expression/model.hpp | 4 +-- test/expression/model_unittest.cpp | 47 ++++++++++++++-------------- 2 files changed, 26 insertions(+), 25 deletions(-) diff --git a/include/autoppl/expression/model.hpp b/include/autoppl/expression/model.hpp index dabb7fc2..264744ac 100644 --- a/include/autoppl/expression/model.hpp +++ b/include/autoppl/expression/model.hpp @@ -32,8 +32,8 @@ struct EqNode { tag_ = tag_cref_; } /* - * Compute pdf of dist_ with value_. - * Assumes that value_ has been assigned with a proper value. + * Compute pdf of underlying distribution with underlying value. + * Assumes that underlying value has been assigned properly. */ dist_value_t pdf() const { return dist_.pdf(tag_.get_value()); } diff --git a/test/expression/model_unittest.cpp b/test/expression/model_unittest.cpp index 62b2c05c..08a81eb8 100644 --- a/test/expression/model_unittest.cpp +++ b/test/expression/model_unittest.cpp @@ -150,35 +150,36 @@ TEST_F(many_tag_dist_fixture, two_tags) y |= MockDist() ); - auto reconfigure = [&](double xv, double yv) { - x.set_value(xv); - y.set_value(yv); - model.update(); - }; - - // both within range xv = 0.2; yv = 1.8; - reconfigure(xv, yv); - EXPECT_EQ(model.pdf(), xv * yv); - EXPECT_EQ(model.log_pdf(), std::log(xv) + std::log(yv)); - // first out of range - xv = 0.000542431; yv = 0.99999; - reconfigure(xv, yv); - EXPECT_EQ(model.pdf(), xv * yv); - EXPECT_EQ(model.log_pdf(), std::log(xv) + std::log(yv)); + x.set_value(xv); + y.set_value(yv); + model.update(); - // second out of range - xv = 0.5142; yv = 2.0000123; - reconfigure(xv, yv); EXPECT_EQ(model.pdf(), xv * yv); EXPECT_EQ(model.log_pdf(), std::log(xv) + std::log(yv)); +} - // both out of range - xv = 14.3; yv = 69.; - reconfigure(xv, yv); - EXPECT_EQ(model.pdf(), xv * yv); - EXPECT_EQ(model.log_pdf(), std::log(xv) + std::log(yv)); +TEST_F(many_tag_dist_fixture, four_tags) +{ + auto model = ( + x |= MockDist(), + y |= MockDist(), + z |= MockDist(), + w |= MockDist() + ); + + xv = 0.2; yv = 1.8; zv = 3.2; wv = 0.3; + + x.set_value(xv); + y.set_value(yv); + z.set_value(zv); + w.set_value(wv); + model.update(); + + EXPECT_EQ(model.pdf(), xv * yv * zv * wv); + EXPECT_EQ(model.log_pdf(), std::log(xv) + std::log(yv) + + std::log(zv) + std::log(wv)); } } // namespace ppl From af1fc7bd920e634c89fdb4d5c78086abc90f82df Mon Sep 17 00:00:00 2001 From: James Yang Date: Wed, 15 Apr 2020 01:36:36 -0400 Subject: [PATCH 09/35] Fix comments about log_pdf --- include/autoppl/expression/model.hpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/include/autoppl/expression/model.hpp b/include/autoppl/expression/model.hpp index 264744ac..50a01195 100644 --- a/include/autoppl/expression/model.hpp +++ b/include/autoppl/expression/model.hpp @@ -39,8 +39,8 @@ struct EqNode { return dist_.pdf(tag_.get_value()); } /* - * Compute log-pdf of dist_ with value x. - * Assumes that value_ has been assigned with a proper value. + * Compute log-pdf of underlying distribution with underlying value. + * Assumes that underlying value has been assigned properly. */ dist_value_t log_pdf() const { return dist_.log_pdf(tag_.get_value()); } From 6ba200c2509d6f83cef873170ba5f25af7a3b391 Mon Sep 17 00:00:00 2001 From: James Yang Date: Wed, 15 Apr 2020 16:54:40 -0400 Subject: [PATCH 10/35] Refactor data from the model expressions to ease mcmc --- include/autoppl/expression/model.hpp | 64 +++++++++++++++++++--------- test/expression/model_unittest.cpp | 45 ++++++++++++++++--- 2 files changed, 85 insertions(+), 24 deletions(-) diff --git a/include/autoppl/expression/model.hpp b/include/autoppl/expression/model.hpp index 50a01195..8007981c 100644 --- a/include/autoppl/expression/model.hpp +++ b/include/autoppl/expression/model.hpp @@ -1,9 +1,21 @@ #pragma once #include #include +#include #include namespace ppl { +namespace details { + +template +struct IdentityTagFunctor +{ + using value_t = typename std::iterator_traits::value_type; + value_t& operator()(value_t& tag) + { return tag; } +}; + +} // namespace details /* * This class represents a "node" in the model expression @@ -18,40 +30,51 @@ struct EqNode EqNode(const tag_t& tag, const dist_t& dist) noexcept - : tag_{} - , tag_cref_{tag} + : orig_tag_cref_{tag} + , comp_tag_cref_{} , dist_{dist} {} /* - * Updates underlying tag by copying from referenced tag. - * This must be called before the model is used, if - * any members of the referenced tag changed. + * Binds computation-required data with this model. + * Underlying reference to computation data for this random variable + * will reference the next data (*begin). + * The next data (*begin) will be initialized with original tag. + * The functor getter MUST return an lvalue reference to the tag. + * TODO: set up compile-time checks for this ^. + * The tag must be the same type as tag_t. */ - void update() - { tag_ = tag_cref_; } + template > + Iter bind_comp_data(Iter begin, Iter end, F getter = F()) + { + assert(begin != end); // MUST have a value to get + getter(*begin) = orig_tag_cref_.get(); // initialize comp data + comp_tag_cref_ = getter(*begin); // set reference to comp data + return ++begin; + } /* * Compute pdf of underlying distribution with underlying value. * Assumes that underlying value has been assigned properly. */ dist_value_t pdf() const - { return dist_.pdf(tag_.get_value()); } + { return dist_.pdf(comp_tag_cref_->get().get_value()); } /* * Compute log-pdf of underlying distribution with underlying value. * Assumes that underlying value has been assigned properly. */ dist_value_t log_pdf() const - { return dist_.log_pdf(tag_.get_value()); } + { return dist_.log_pdf(comp_tag_cref_->get().get_value()); } private: using tag_cref_t = std::reference_wrapper; + using opt_tag_cref_t = std::optional; - tag_t tag_; // cache optimization - tag_cref_t tag_cref_; // (const) reference of the tag since any configuration - // may be changed until right before update - dist_t dist_; // distribution associated with tag + tag_cref_t orig_tag_cref_; // (const) reference of the original tag since + // any configuration may be changed until right before update + opt_tag_cref_t comp_tag_cref_; // reference the tag needed in computation + dist_t dist_; // distribution associated with tag }; /* @@ -69,18 +92,21 @@ struct GlueNode >; GlueNode(const left_node_t& lhs, - const right_node_t& rhs) + const right_node_t& rhs) noexcept : left_node_{lhs} , right_node_{rhs} {} /* - * Updates left node first then right node by calling "update" on each. - * This must be called before the model is used, - * if either left or right must be updated. + * Binds computational data in order from left to right. + * In other words, same order as user would list the model expressions. */ - void update() - { left_node_.update(); right_node_.update(); } + template > + Iter bind_comp_data(Iter begin, Iter end, F getter = F()) + { + Iter new_begin = left_node_.bind_comp_data(begin, end, getter); + return right_node_.bind_comp_data(new_begin, end, getter); + } /* * Computes left node joint pdf then right node joint pdf diff --git a/test/expression/model_unittest.cpp b/test/expression/model_unittest.cpp index 08a81eb8..551675a2 100644 --- a/test/expression/model_unittest.cpp +++ b/test/expression/model_unittest.cpp @@ -44,12 +44,12 @@ struct MockDist }; /* - * Fixture for testing one rv_tag with distribution. + * Fixture for testing one tag with distribution. */ struct tag_dist_fixture : ::testing::Test { protected: - MockTag x; + MockTag x, comp_data; using model_t = std::decay_t; model_t model = (x |= MockDist()); double val; @@ -57,7 +57,8 @@ struct tag_dist_fixture : ::testing::Test void reconfigure(double val) { x.set_value(val); - model.update(); + auto ptr = model.bind_comp_data(&comp_data, &comp_data + 1); + EXPECT_EQ(ptr, &comp_data + 1); } }; @@ -136,11 +137,19 @@ TEST_F(tag_dist_fixture, log_pdf_invalid) // Model with many RV (no dependencies) TESTS ////////////////////////////////////////////////////// +/* + * Fixture for testing many tags with distributions. + */ struct many_tag_dist_fixture : ::testing::Test { protected: + std::vector comp_data; MockTag x, y, z, w; double xv, yv, zv, wv; + + many_tag_dist_fixture() + : comp_data(4) + {} }; TEST_F(many_tag_dist_fixture, two_tags) @@ -154,7 +163,7 @@ TEST_F(many_tag_dist_fixture, two_tags) x.set_value(xv); y.set_value(yv); - model.update(); + model.bind_comp_data(comp_data.begin(), std::next(comp_data.begin(), 2)); EXPECT_EQ(model.pdf(), xv * yv); EXPECT_EQ(model.log_pdf(), std::log(xv) + std::log(yv)); @@ -175,11 +184,37 @@ TEST_F(many_tag_dist_fixture, four_tags) y.set_value(yv); z.set_value(zv); w.set_value(wv); - model.update(); + model.bind_comp_data(comp_data.begin(), comp_data.end()); EXPECT_EQ(model.pdf(), xv * yv * zv * wv); EXPECT_EQ(model.log_pdf(), std::log(xv) + std::log(yv) + std::log(zv) + std::log(wv)); } +TEST_F(many_tag_dist_fixture, four_tags_correct_bind) +{ + auto model = ( + x |= MockDist(), + y |= MockDist(), + z |= MockDist(), + w |= MockDist() + ); + + xv = 0.2; yv = 1.8; zv = 3.2; wv = 0.3; + + x.set_value(xv); + y.set_value(yv); + z.set_value(zv); + w.set_value(wv); + model.bind_comp_data(comp_data.begin(), comp_data.end()); + + // test that the computation data was initialized in the same + // order as the variables listed in the model. + + EXPECT_EQ(comp_data[0].get_value(), xv); + EXPECT_EQ(comp_data[1].get_value(), yv); + EXPECT_EQ(comp_data[2].get_value(), zv); + EXPECT_EQ(comp_data[3].get_value(), wv); +} + } // namespace ppl From d0158b0d8705a1fa531825dfed3898a19c206122 Mon Sep 17 00:00:00 2001 From: James Yang Date: Fri, 17 Apr 2020 15:11:09 -0400 Subject: [PATCH 11/35] Add design document --- doc/design/model_inttest.cpp | 106 +++++++++++++++++++++++++++++++++++ 1 file changed, 106 insertions(+) create mode 100644 doc/design/model_inttest.cpp diff --git a/doc/design/model_inttest.cpp b/doc/design/model_inttest.cpp new file mode 100644 index 00000000..64caf2ca --- /dev/null +++ b/doc/design/model_inttest.cpp @@ -0,0 +1,106 @@ +#include "gtest/gtest.h" +#include +#include +#include + +namespace ppl { + +template +struct BracketNode +{ + VectorType v; + IndexType i; +}; + +struct myvector +{ + + rv_tag operator[](rv_tag) + { + return rv + } + std::vector v; // 3 things +}; + +template +auto normal(const MuType& mu, const SigType& sig) +{ + Normal(mu, sig); +} + +TEST(dummy, dummy_test) +{ + double x_data = 2.3; // 1-sample data + + std::vector sampled_theta_1(100); + std::vector sampled_theta_2(100); + + double* ptr; + rv_tag x; + rv_tag theta_1(sampled_theta_1.data()); + rv_tag theta_2(sampled_theta_2.data()); + + std::vector> v; + std::for_each(..., ... , [](){v[i].set_sample_storage(&mat.row(i));}); + + x.observe(x_data); + + x_1.observe(...); + x_2.observe(...); + + auto model = ( + mu |= uniform(-10000, 10000), + y |= uniform({1,2,3}) // + x_1 |= normal(mu[y], 1), + x_2 |= normal(mu[y], 1), + ); + + x.observe(...); + + rv_tag var, mu, x; + auto normal_model = ( + var |= normal(0,1), + mu |= normal(1,5), + x |= normal(mu, var) + ); + + std::vector var_storage(1000); + std::vector mu_storage(1000); + + var.set_storage(var_storage.data()); + mu.set_storage(mu_storage.data()); + + metropolis_hastings(model, 1000, 400); + + auto gmm_model = ( + mu |= + ); + + std::vector> vec(model.param_num); + model.bind_storage(vec.begin(), vec.end(), ...); + model.pdf(); + + metropolis_hastings(model, 100); + + std::vector sampled_theta_1_again(1000); + std::vector sampled_theta_2_again(1000); + + theta_1.set_storage(sampled_theta_1_again.data()); + theta_2.set_storage(sampled_theta_2_again.data()); + + metropolis_hastings(model, 1000); + + + + + + + + auto model = ( + w |= normal(0,1), + y |= normal(w*x, 1) + ) + metropolis_hastings(modeli) +} + +} From 60ce086c58a805216a311ee4684ffd660a5c2983 Mon Sep 17 00:00:00 2001 From: James Yang Date: Sat, 18 Apr 2020 12:30:07 -0400 Subject: [PATCH 12/35] Simplify unittest logic --- test/expression/model_unittest.cpp | 14 ++++++++++---- 1 file changed, 10 insertions(+), 4 deletions(-) diff --git a/test/expression/model_unittest.cpp b/test/expression/model_unittest.cpp index 551675a2..ecff9c24 100644 --- a/test/expression/model_unittest.cpp +++ b/test/expression/model_unittest.cpp @@ -1,4 +1,5 @@ #include +#include #include "gtest/gtest.h" #include @@ -49,16 +50,21 @@ struct MockDist struct tag_dist_fixture : ::testing::Test { protected: - MockTag x, comp_data; + MockTag x; + std::array comp_data; using model_t = std::decay_t; model_t model = (x |= MockDist()); double val; + tag_dist_fixture() + { + auto ptr = model.bind_comp_data(comp_data.begin(), comp_data.end()); + EXPECT_EQ(ptr, comp_data.end()); + } + void reconfigure(double val) { - x.set_value(val); - auto ptr = model.bind_comp_data(&comp_data, &comp_data + 1); - EXPECT_EQ(ptr, &comp_data + 1); + comp_data[0].set_value(val); } }; From 5c0a4e331282f62facaee6976da50e491b3d5139 Mon Sep 17 00:00:00 2001 From: James Yang Date: Sat, 18 Apr 2020 12:52:31 -0400 Subject: [PATCH 13/35] Simplify model pdf to read from original tags --- include/autoppl/expression/model.hpp | 35 ++--------------------- test/expression/model_unittest.cpp | 42 +--------------------------- 2 files changed, 3 insertions(+), 74 deletions(-) diff --git a/include/autoppl/expression/model.hpp b/include/autoppl/expression/model.hpp index 8007981c..ae236a4e 100644 --- a/include/autoppl/expression/model.hpp +++ b/include/autoppl/expression/model.hpp @@ -31,41 +31,22 @@ struct EqNode EqNode(const tag_t& tag, const dist_t& dist) noexcept : orig_tag_cref_{tag} - , comp_tag_cref_{} , dist_{dist} {} - /* - * Binds computation-required data with this model. - * Underlying reference to computation data for this random variable - * will reference the next data (*begin). - * The next data (*begin) will be initialized with original tag. - * The functor getter MUST return an lvalue reference to the tag. - * TODO: set up compile-time checks for this ^. - * The tag must be the same type as tag_t. - */ - template > - Iter bind_comp_data(Iter begin, Iter end, F getter = F()) - { - assert(begin != end); // MUST have a value to get - getter(*begin) = orig_tag_cref_.get(); // initialize comp data - comp_tag_cref_ = getter(*begin); // set reference to comp data - return ++begin; - } - /* * Compute pdf of underlying distribution with underlying value. * Assumes that underlying value has been assigned properly. */ dist_value_t pdf() const - { return dist_.pdf(comp_tag_cref_->get().get_value()); } + { return dist_.pdf(orig_tag_cref_.get().get_value()); } /* * Compute log-pdf of underlying distribution with underlying value. * Assumes that underlying value has been assigned properly. */ dist_value_t log_pdf() const - { return dist_.log_pdf(comp_tag_cref_->get().get_value()); } + { return dist_.log_pdf(orig_tag_cref_.get().get_value()); } private: using tag_cref_t = std::reference_wrapper; @@ -73,7 +54,6 @@ struct EqNode tag_cref_t orig_tag_cref_; // (const) reference of the original tag since // any configuration may be changed until right before update - opt_tag_cref_t comp_tag_cref_; // reference the tag needed in computation dist_t dist_; // distribution associated with tag }; @@ -97,17 +77,6 @@ struct GlueNode , right_node_{rhs} {} - /* - * Binds computational data in order from left to right. - * In other words, same order as user would list the model expressions. - */ - template > - Iter bind_comp_data(Iter begin, Iter end, F getter = F()) - { - Iter new_begin = left_node_.bind_comp_data(begin, end, getter); - return right_node_.bind_comp_data(new_begin, end, getter); - } - /* * Computes left node joint pdf then right node joint pdf * and returns the product of the two. diff --git a/test/expression/model_unittest.cpp b/test/expression/model_unittest.cpp index ecff9c24..7be5171f 100644 --- a/test/expression/model_unittest.cpp +++ b/test/expression/model_unittest.cpp @@ -51,20 +51,13 @@ struct tag_dist_fixture : ::testing::Test { protected: MockTag x; - std::array comp_data; using model_t = std::decay_t; model_t model = (x |= MockDist()); double val; - tag_dist_fixture() - { - auto ptr = model.bind_comp_data(comp_data.begin(), comp_data.end()); - EXPECT_EQ(ptr, comp_data.end()); - } - void reconfigure(double val) { - comp_data[0].set_value(val); + x.set_value(val); } }; @@ -149,13 +142,8 @@ TEST_F(tag_dist_fixture, log_pdf_invalid) struct many_tag_dist_fixture : ::testing::Test { protected: - std::vector comp_data; MockTag x, y, z, w; double xv, yv, zv, wv; - - many_tag_dist_fixture() - : comp_data(4) - {} }; TEST_F(many_tag_dist_fixture, two_tags) @@ -169,7 +157,6 @@ TEST_F(many_tag_dist_fixture, two_tags) x.set_value(xv); y.set_value(yv); - model.bind_comp_data(comp_data.begin(), std::next(comp_data.begin(), 2)); EXPECT_EQ(model.pdf(), xv * yv); EXPECT_EQ(model.log_pdf(), std::log(xv) + std::log(yv)); @@ -190,37 +177,10 @@ TEST_F(many_tag_dist_fixture, four_tags) y.set_value(yv); z.set_value(zv); w.set_value(wv); - model.bind_comp_data(comp_data.begin(), comp_data.end()); EXPECT_EQ(model.pdf(), xv * yv * zv * wv); EXPECT_EQ(model.log_pdf(), std::log(xv) + std::log(yv) + std::log(zv) + std::log(wv)); } -TEST_F(many_tag_dist_fixture, four_tags_correct_bind) -{ - auto model = ( - x |= MockDist(), - y |= MockDist(), - z |= MockDist(), - w |= MockDist() - ); - - xv = 0.2; yv = 1.8; zv = 3.2; wv = 0.3; - - x.set_value(xv); - y.set_value(yv); - z.set_value(zv); - w.set_value(wv); - model.bind_comp_data(comp_data.begin(), comp_data.end()); - - // test that the computation data was initialized in the same - // order as the variables listed in the model. - - EXPECT_EQ(comp_data[0].get_value(), xv); - EXPECT_EQ(comp_data[1].get_value(), yv); - EXPECT_EQ(comp_data[2].get_value(), zv); - EXPECT_EQ(comp_data[3].get_value(), wv); -} - } // namespace ppl From d4ce9ae1d83d40b9f3c6f63488312e3c48ac1885 Mon Sep 17 00:00:00 2001 From: Jacob Austin Date: Sat, 18 Apr 2020 13:11:43 -0400 Subject: [PATCH 14/35] updated uniform to support rv_tag parameters --- doc/design/model_design2.cpp | 32 ++++++++++++++++++++++++++ include/autoppl/expression/rv_tag.hpp | 5 ++++ include/autoppl/expression/uniform.hpp | 32 ++++++++++++++++++-------- test/CMakeLists.txt | 1 + test/expression/uniform_unittest.cpp | 29 +++++++++++++++++++++++ 5 files changed, 90 insertions(+), 9 deletions(-) create mode 100644 doc/design/model_design2.cpp create mode 100644 test/expression/uniform_unittest.cpp diff --git a/doc/design/model_design2.cpp b/doc/design/model_design2.cpp new file mode 100644 index 00000000..5aad3057 --- /dev/null +++ b/doc/design/model_design2.cpp @@ -0,0 +1,32 @@ +Y ~ W.x + epsilon +Y ~ N(W.x, sigma^2) + +Parameter X {4.0}; // observed +Parameter Y {5.0}; // observed +Parameter W; // hidden +​ +Model m1 = Model( // Model class defines a distribution over existing Parameters. + W |= Uniform(-10, 10), // linear regression + Y |= Normal(W * X, 3), // overload multiplication to build a graph from W * X +​); + +Model m2 = Model( + W |= Normal(0, 1), // ridge regression instead + Y |= Normal(W * X, 3), +); +​ +m1.sample(1000); + +(3*x).pdf(10) => x.pdf(10 / 3) + +X.observe(3); // observe more data + +// P(Y, W | X) = P(Y | W, X) P(W | X) which is doable for multiple samples, just need to +// assert len(Y) == len(X) and then multiply out over all pairs of (X, Y) values. + +// P(Y | X) => this is a fine distribution, but I can't talk about P(Y, X) or P(X | Y) until I put a prior on Y. +// I don't have a joint distribution yet. + +// Some issues: +// how do we do (x ** 2).pdf(5)? This is pretty damn hard for non-bijective functions, need to integrate? +// \ No newline at end of file diff --git a/include/autoppl/expression/rv_tag.hpp b/include/autoppl/expression/rv_tag.hpp index 00ccaa65..19063603 100644 --- a/include/autoppl/expression/rv_tag.hpp +++ b/include/autoppl/expression/rv_tag.hpp @@ -37,6 +37,9 @@ struct rv_tag : rv_tag(0, storage_ptr) {} + rv_tag(value_t value) noexcept + : rv_tag(value, nullptr) {} + rv_tag() noexcept : rv_tag(0, nullptr) {} @@ -51,6 +54,8 @@ struct rv_tag void set_state(state_t state) { state_ = state; } state_t get_state() const { return state_; } + operator value_t () const { return value_; } + /* * Sets underlying value to "value". * Additionally modifies the tag to be considered as data. diff --git a/include/autoppl/expression/uniform.hpp b/include/autoppl/expression/uniform.hpp index 562d2c6b..15547ace 100644 --- a/include/autoppl/expression/uniform.hpp +++ b/include/autoppl/expression/uniform.hpp @@ -8,37 +8,51 @@ namespace ppl { // TODO: change name to UniformDist and make class template. // uniform should be a function that creates this kind of object. -struct uniform + +template +struct Uniform { using value_t = double; using dist_value_t = double; - uniform(value_t min, value_t max) - : min_{min}, max_{max} - { assert(min_ < max_); } + Uniform(min_type min, max_type max) + : min_{min}, max_{max} { assert(static_cast(min_) < static_cast(max_)); } // TODO: tag this class as "TriviallySamplable"? template value_t sample(GeneratorType& gen) const { - std::uniform_real_distribution dist(min_, max_); + value_t min, max; + min = static_cast(min_); + max = static_cast(max_); + + std::uniform_real_distribution dist(min, max); return dist(gen); } dist_value_t pdf(value_t x) const { - return (min_ < x && x < max_) ? 1./(max_ - min_) : 0; + value_t min, max; + min = static_cast(min_); + max = static_cast(max_); + + return (min < x && x < max) ? 1. / (max - min) : 0; } dist_value_t log_pdf(value_t x) const { - return (min_ < x && x < max_) ? - -std::log(max_ - min_) : + value_t min, max; + min = static_cast(min_); + max = static_cast(max_); + + return (min < x && x < max) ? + -std::log(max - min) : std::numeric_limits::lowest(); } private: - value_t min_, max_; + min_type min_; + max_type max_; }; } // namespace ppl diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index ef3cd0e2..be61815b 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -11,6 +11,7 @@ endif() add_executable(expression_unittest ${CMAKE_CURRENT_SOURCE_DIR}/expression/model_unittest.cpp + ${CMAKE_CURRENT_SOURCE_DIR}/expression/uniform_unittest.cpp ) target_compile_options(expression_unittest PRIVATE -g -Wall -Werror -Wextra) target_include_directories(expression_unittest PRIVATE ${GTEST_DIR}/include) diff --git a/test/expression/uniform_unittest.cpp b/test/expression/uniform_unittest.cpp new file mode 100644 index 00000000..99b9448c --- /dev/null +++ b/test/expression/uniform_unittest.cpp @@ -0,0 +1,29 @@ +#include +#include +#include + +#include +#include + +#include "gtest/gtest.h" + +namespace ppl { + +struct uniform_dist_fixture : ::testing::Test { +protected: + rv_tag x {0.5}; + rv_tag y {0.1}; + Uniform dist1 = Uniform(0., 1.); + Uniform > dist2 = Uniform(0., x); + Uniform, rv_tag > dist3 = Uniform(y, x); +}; + +TEST_F(uniform_dist_fixture, simple_uniform) { + EXPECT_EQ(dist1.pdf(1.1), 0.0); + EXPECT_EQ(dist2.pdf(1.0), 0.0); + EXPECT_EQ(dist2.pdf(0.25), 2.0); + EXPECT_EQ(dist3.pdf(-0.1), 0.0); + EXPECT_EQ(dist3.pdf(0.25), 2.5); +} + +} // ppl \ No newline at end of file From d8170ea5d9b14feffa9b8d155ad29e7f4211d7c1 Mon Sep 17 00:00:00 2001 From: Jacob Austin Date: Sat, 18 Apr 2020 13:35:38 -0400 Subject: [PATCH 15/35] renamed var --- include/autoppl/expression/model.hpp | 36 +++++++++---------- include/autoppl/expression/traits.hpp | 13 ++++--- .../expression/{rv_tag.hpp => variable.hpp} | 36 +++++++++---------- test/CMakeLists.txt | 16 ++++++++- .../uniform_unittest.cpp | 13 +++---- test/expression/model_unittest.cpp | 30 ++++++++-------- 6 files changed, 81 insertions(+), 63 deletions(-) rename include/autoppl/expression/{rv_tag.hpp => variable.hpp} (66%) rename test/{expression => distribution}/uniform_unittest.cpp (64%) diff --git a/include/autoppl/expression/model.hpp b/include/autoppl/expression/model.hpp index ae236a4e..a9566e8f 100644 --- a/include/autoppl/expression/model.hpp +++ b/include/autoppl/expression/model.hpp @@ -8,29 +8,29 @@ namespace ppl { namespace details { template -struct IdentityTagFunctor +struct IdentityVarFunctor { using value_t = typename std::iterator_traits::value_type; - value_t& operator()(value_t& tag) - { return tag; } + value_t& operator()(value_t& var) + { return var; } }; } // namespace details /* * This class represents a "node" in the model expression - * that relates a tag with a distribution. + * that relates a var with a distribution. */ -template +template struct EqNode { - using tag_t = TagType; + using var_t = VarType; using dist_t = DistType; using dist_value_t = typename dist_traits::dist_value_t; - EqNode(const tag_t& tag, + EqNode(const var_t& var, const dist_t& dist) noexcept - : orig_tag_cref_{tag} + : orig_var_cref_{var} , dist_{dist} {} @@ -39,22 +39,22 @@ struct EqNode * Assumes that underlying value has been assigned properly. */ dist_value_t pdf() const - { return dist_.pdf(orig_tag_cref_.get().get_value()); } + { return dist_.pdf(orig_var_cref_.get().get_value()); } /* * Compute log-pdf of underlying distribution with underlying value. * Assumes that underlying value has been assigned properly. */ dist_value_t log_pdf() const - { return dist_.log_pdf(orig_tag_cref_.get().get_value()); } + { return dist_.log_pdf(orig_var_cref_.get().get_value()); } private: - using tag_cref_t = std::reference_wrapper; - using opt_tag_cref_t = std::optional; + using var_cref_t = std::reference_wrapper; + using opt_var_cref_t = std::optional; - tag_cref_t orig_tag_cref_; // (const) reference of the original tag since + var_cref_t orig_var_cref_; // (const) reference of the original var since // any configuration may be changed until right before update - dist_t dist_; // distribution associated with tag + dist_t dist_; // distribution associated with var }; /* @@ -104,14 +104,14 @@ struct GlueNode // with concepts! /* - * Builds an EqNode to associate tag with dist. + * Builds an EqNode to associate var with dist. * Ex. x |= uniform(0,1) */ -template -constexpr inline auto operator|=(const TagType& tag, +template +constexpr inline auto operator|=(const VarType& var, const DistType& dist) { - return EqNode(tag, dist); + return EqNode(var, dist); } /* diff --git a/include/autoppl/expression/traits.hpp b/include/autoppl/expression/traits.hpp index dce5efc6..bb4c83d2 100644 --- a/include/autoppl/expression/traits.hpp +++ b/include/autoppl/expression/traits.hpp @@ -1,4 +1,5 @@ #pragma once +#include namespace ppl { @@ -8,12 +9,14 @@ namespace ppl { * Users should rely on these classes to grab member aliases. */ -template -struct tag_traits +template +struct var_traits { - using value_t = typename TagType::value_t; - using pointer_t = typename TagType::pointer_t; - using state_t = typename TagType::state_t; + using value_t = typename VarType::value_t; + using pointer_t = typename VarType::pointer_t; + using state_t = typename VarType::state_t; + + static_assert(std::is_convertible_v); }; template diff --git a/include/autoppl/expression/rv_tag.hpp b/include/autoppl/expression/variable.hpp similarity index 66% rename from include/autoppl/expression/rv_tag.hpp rename to include/autoppl/expression/variable.hpp index 19063603..3c332e9f 100644 --- a/include/autoppl/expression/rv_tag.hpp +++ b/include/autoppl/expression/variable.hpp @@ -3,45 +3,45 @@ namespace ppl { /* - * The possible states for a tag. - * By default, all tags should be considered as a parameter. + * The possible states for a var. + * By default, all vars should be considered as a parameter. * TODO: maybe move in a different file? */ -enum class tag_state : bool { +enum class var_state : bool { data, parameter }; /* - * rv_tag is a light-weight structure that represents a univariate random variable. + * Variable is a light-weight structure that represents a univariate random variable. * It acts as an intermediate layer of communication between - * a model expression and the users, who must supply storage of values associated with this tag. + * a model expression and the users, who must supply storage of values associated with this var. */ template -struct rv_tag +struct Variable { using value_t = ValueType; using pointer_t = value_t*; using const_pointer_t = const value_t*; - using state_t = tag_state; + using state_t = var_state; // constructors - rv_tag(value_t value, + Variable(value_t value, pointer_t storage_ptr) noexcept : value_{value} , storage_ptr_{storage_ptr} , state_{state_t::parameter} {} - rv_tag(pointer_t storage_ptr) noexcept - : rv_tag(0, storage_ptr) + Variable(pointer_t storage_ptr) noexcept + : Variable(0, storage_ptr) {} - rv_tag(value_t value) noexcept - : rv_tag(value, nullptr) {} + Variable(value_t value) noexcept + : Variable(value, nullptr) {} - rv_tag() noexcept - : rv_tag(0, nullptr) + Variable() noexcept + : Variable(0, nullptr) {} void set_value(value_t value) { value_ = value; } @@ -58,7 +58,7 @@ struct rv_tag /* * Sets underlying value to "value". - * Additionally modifies the tag to be considered as data. + * Additionally modifies the var to be considered as data. * Equivalent to calling set_value(value) then set_state(state). */ void observe(value_t value) @@ -68,14 +68,14 @@ struct rv_tag } private: - value_t value_; // store value associated with tag + value_t value_; // store value associated with var pointer_t storage_ptr_; // points to beginning of storage // storage is assumed to be contiguous state_t state_; // state to determine if data or param }; // Useful aliases -using cont_rv_tag = rv_tag; // continuous RV tag -using disc_rv_tag = rv_tag; // discrete RV tag +using cont_var = Variable; // continuous RV var +using disc_var = Variable; // discrete RV var } // namespace ppl diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index be61815b..998e9576 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -11,7 +11,6 @@ endif() add_executable(expression_unittest ${CMAKE_CURRENT_SOURCE_DIR}/expression/model_unittest.cpp - ${CMAKE_CURRENT_SOURCE_DIR}/expression/uniform_unittest.cpp ) target_compile_options(expression_unittest PRIVATE -g -Wall -Werror -Wextra) target_include_directories(expression_unittest PRIVATE ${GTEST_DIR}/include) @@ -20,3 +19,18 @@ if (AUTOPPL_ENABLE_TEST_COVERAGE) endif() target_link_libraries(expression_unittest gtest_main pthread ${PROJECT_NAME}) add_test(expression_unittest expression_unittest) + +###################################################### +# Distribution Test +###################################################### + +add_executable(distribution_unittest + ${CMAKE_CURRENT_SOURCE_DIR}/distribution/uniform_unittest.cpp + ) +target_compile_options(distribution_unittest PRIVATE -g -Wall -Werror -Wextra) +target_include_directories(distribution_unittest PRIVATE ${GTEST_DIR}/include) +if (AUTOPPL_ENABLE_TEST_COVERAGE) + target_link_libraries(distribution_unittest gcov) +endif() +target_link_libraries(distribution_unittest gtest_main pthread ${PROJECT_NAME}) +add_test(distribution_unittest distribution_unittest) \ No newline at end of file diff --git a/test/expression/uniform_unittest.cpp b/test/distribution/uniform_unittest.cpp similarity index 64% rename from test/expression/uniform_unittest.cpp rename to test/distribution/uniform_unittest.cpp index 99b9448c..f9f52a26 100644 --- a/test/expression/uniform_unittest.cpp +++ b/test/distribution/uniform_unittest.cpp @@ -1,6 +1,5 @@ #include -#include -#include +#include #include #include @@ -11,17 +10,19 @@ namespace ppl { struct uniform_dist_fixture : ::testing::Test { protected: - rv_tag x {0.5}; - rv_tag y {0.1}; + Variable x {0.5}; + Variable y {0.1}; Uniform dist1 = Uniform(0., 1.); - Uniform > dist2 = Uniform(0., x); - Uniform, rv_tag > dist3 = Uniform(y, x); + Uniform > dist2 = Uniform(0., x); + Uniform, Variable > dist3 = Uniform(y, x); }; TEST_F(uniform_dist_fixture, simple_uniform) { EXPECT_EQ(dist1.pdf(1.1), 0.0); + EXPECT_EQ(dist2.pdf(1.0), 0.0); EXPECT_EQ(dist2.pdf(0.25), 2.0); + EXPECT_EQ(dist3.pdf(-0.1), 0.0); EXPECT_EQ(dist3.pdf(0.25), 2.5); } diff --git a/test/expression/model_unittest.cpp b/test/expression/model_unittest.cpp index 7be5171f..9007b3ab 100644 --- a/test/expression/model_unittest.cpp +++ b/test/expression/model_unittest.cpp @@ -10,10 +10,10 @@ namespace ppl { ////////////////////////////////////////////////////// /* - * Mock tag object for testing purposes. - * Must meet some of the requirements of actual tag types. + * Mock var object for testing purposes. + * Must meet some of the requirements of actual var types. */ -struct MockTag +struct MockVar { using value_t = double; using pointer_t = double*; @@ -45,12 +45,12 @@ struct MockDist }; /* - * Fixture for testing one tag with distribution. + * Fixture for testing one var with distribution. */ -struct tag_dist_fixture : ::testing::Test +struct var_dist_fixture : ::testing::Test { protected: - MockTag x; + MockVar x; using model_t = std::decay_t; model_t model = (x |= MockDist()); double val; @@ -61,7 +61,7 @@ struct tag_dist_fixture : ::testing::Test } }; -TEST_F(tag_dist_fixture, pdf_valid) +TEST_F(var_dist_fixture, pdf_valid) { // MockDist pdf is identity function // so we may simply compare model.pdf() with val. @@ -79,7 +79,7 @@ TEST_F(tag_dist_fixture, pdf_valid) EXPECT_EQ(model.pdf(), val); } -TEST_F(tag_dist_fixture, pdf_invalid) +TEST_F(var_dist_fixture, pdf_invalid) { val = 0.000004123; reconfigure(val); @@ -98,7 +98,7 @@ TEST_F(tag_dist_fixture, pdf_invalid) EXPECT_EQ(model.pdf(), val); } -TEST_F(tag_dist_fixture, log_pdf_valid) +TEST_F(var_dist_fixture, log_pdf_valid) { val = 0.000001; reconfigure(val); @@ -113,7 +113,7 @@ TEST_F(tag_dist_fixture, log_pdf_valid) EXPECT_EQ(model.log_pdf(), std::log(val)); } -TEST_F(tag_dist_fixture, log_pdf_invalid) +TEST_F(var_dist_fixture, log_pdf_invalid) { val = 0.000004123; reconfigure(val); @@ -137,16 +137,16 @@ TEST_F(tag_dist_fixture, log_pdf_invalid) ////////////////////////////////////////////////////// /* - * Fixture for testing many tags with distributions. + * Fixture for testing many vars with distributions. */ -struct many_tag_dist_fixture : ::testing::Test +struct many_var_dist_fixture : ::testing::Test { protected: - MockTag x, y, z, w; + MockVar x, y, z, w; double xv, yv, zv, wv; }; -TEST_F(many_tag_dist_fixture, two_tags) +TEST_F(many_var_dist_fixture, two_vars) { auto model = ( x |= MockDist(), @@ -162,7 +162,7 @@ TEST_F(many_tag_dist_fixture, two_tags) EXPECT_EQ(model.log_pdf(), std::log(xv) + std::log(yv)); } -TEST_F(many_tag_dist_fixture, four_tags) +TEST_F(many_var_dist_fixture, four_vars) { auto model = ( x |= MockDist(), From a9f0a04a612cc07cd46b82412e8ab45f57de6c3e Mon Sep 17 00:00:00 2001 From: Jacob Austin Date: Sat, 18 Apr 2020 13:38:58 -0400 Subject: [PATCH 16/35] moved distributions to a separate folder --- include/autoppl/{expression => distribution}/uniform.hpp | 0 test/distribution/uniform_unittest.cpp | 2 +- 2 files changed, 1 insertion(+), 1 deletion(-) rename include/autoppl/{expression => distribution}/uniform.hpp (100%) diff --git a/include/autoppl/expression/uniform.hpp b/include/autoppl/distribution/uniform.hpp similarity index 100% rename from include/autoppl/expression/uniform.hpp rename to include/autoppl/distribution/uniform.hpp diff --git a/test/distribution/uniform_unittest.cpp b/test/distribution/uniform_unittest.cpp index f9f52a26..f9b98e99 100644 --- a/test/distribution/uniform_unittest.cpp +++ b/test/distribution/uniform_unittest.cpp @@ -1,4 +1,4 @@ -#include +#include #include #include From a371f33fbd3408213cfc8451684d487cfc43eb3a Mon Sep 17 00:00:00 2001 From: Jacob Austin Date: Sat, 18 Apr 2020 14:09:26 -0400 Subject: [PATCH 17/35] added normal distributions --- include/autoppl/distribution/normal.hpp | 58 +++++++++++++++++++++++++ test/CMakeLists.txt | 3 +- test/distribution/normal_unittest.cpp | 29 +++++++++++++ test/distribution/uniform_unittest.cpp | 10 ++--- 4 files changed, 94 insertions(+), 6 deletions(-) create mode 100644 include/autoppl/distribution/normal.hpp create mode 100644 test/distribution/normal_unittest.cpp diff --git a/include/autoppl/distribution/normal.hpp b/include/autoppl/distribution/normal.hpp new file mode 100644 index 00000000..8e1e6181 --- /dev/null +++ b/include/autoppl/distribution/normal.hpp @@ -0,0 +1,58 @@ +#pragma once +#include + +#include +#include +#include +#include + +namespace ppl { + +// TODO: change name to NormalDist and make class template. +// normal should be a function that creates this kind of object. + +template +struct Normal { + using value_t = double; + using dist_value_t = double; + + static_assert(std::is_convertible_v); + static_assert(std::is_convertible_v); + + Normal(mean_type mean, var_type var) + : mean_{mean}, var_{var} { + assert(static_cast(var_) > 0); + }; + + template + value_t sample(GeneratorType& gen) const { + value_t mean, var; + mean = static_cast(mean_); + var = static_cast(var_); + + std::normal_distribution dist(mean, var); + return dist(gen); + } + + dist_value_t pdf(value_t x) const { + value_t mean, var; + mean = static_cast(mean_); + var = static_cast(var_); + + return std::exp(- 0.5 * std::pow(x - mean, 2) / var) / (std::sqrt(var * 2 * M_PI)); + } + + dist_value_t log_pdf(value_t x) const { + value_t mean, var; + mean = static_cast(mean_); + var = static_cast(var_); + + return (-0.5 * std::pow(x - mean, 2) / var) - 0.5 * (std::log(var) + std::log(2) + std::log(M_PI)); + } + + private: + mean_type mean_; + var_type var_; +}; + +} // namespace ppl diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 998e9576..7f93364b 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -26,6 +26,7 @@ add_test(expression_unittest expression_unittest) add_executable(distribution_unittest ${CMAKE_CURRENT_SOURCE_DIR}/distribution/uniform_unittest.cpp + ${CMAKE_CURRENT_SOURCE_DIR}/distribution/normal_unittest.cpp ) target_compile_options(distribution_unittest PRIVATE -g -Wall -Werror -Wextra) target_include_directories(distribution_unittest PRIVATE ${GTEST_DIR}/include) @@ -33,4 +34,4 @@ if (AUTOPPL_ENABLE_TEST_COVERAGE) target_link_libraries(distribution_unittest gcov) endif() target_link_libraries(distribution_unittest gtest_main pthread ${PROJECT_NAME}) -add_test(distribution_unittest distribution_unittest) \ No newline at end of file +add_test(distribution_unittest distribution_unittest) diff --git a/test/distribution/normal_unittest.cpp b/test/distribution/normal_unittest.cpp new file mode 100644 index 00000000..7beb42b3 --- /dev/null +++ b/test/distribution/normal_unittest.cpp @@ -0,0 +1,29 @@ +#include +#include + +#include +#include + +#include "gtest/gtest.h" + +namespace ppl { + +struct normal_dist_fixture : ::testing::Test { +protected: + Variable mu {0.}; + Variable sigma {1.}; + Normal dist1 = Normal(0., 1.); + Normal, Variable > dist2 = Normal(mu, sigma); +}; + +TEST_F(normal_dist_fixture, simple_gaussian) { + EXPECT_DOUBLE_EQ(dist1.pdf(0.0), 0.3989422804014327); + EXPECT_DOUBLE_EQ(dist1.pdf(-0.5), 0.3520653267642995); + EXPECT_DOUBLE_EQ(dist1.pdf(4), 0.00013383022576488537); + + EXPECT_DOUBLE_EQ(dist1.log_pdf(0.0), std::log(dist1.pdf(0.0))); + EXPECT_DOUBLE_EQ(dist1.log_pdf(-0.5), std::log(dist1.pdf(-0.5))); + EXPECT_DOUBLE_EQ(dist1.log_pdf(4), std::log(dist1.pdf(4))); +} + +} // ppl \ No newline at end of file diff --git a/test/distribution/uniform_unittest.cpp b/test/distribution/uniform_unittest.cpp index f9b98e99..fa13c601 100644 --- a/test/distribution/uniform_unittest.cpp +++ b/test/distribution/uniform_unittest.cpp @@ -18,13 +18,13 @@ struct uniform_dist_fixture : ::testing::Test { }; TEST_F(uniform_dist_fixture, simple_uniform) { - EXPECT_EQ(dist1.pdf(1.1), 0.0); + EXPECT_DOUBLE_EQ(dist1.pdf(1.1), 0.0); - EXPECT_EQ(dist2.pdf(1.0), 0.0); - EXPECT_EQ(dist2.pdf(0.25), 2.0); + EXPECT_DOUBLE_EQ(dist2.pdf(1.0), 0.0); + EXPECT_DOUBLE_EQ(dist2.pdf(0.25), 2.0); - EXPECT_EQ(dist3.pdf(-0.1), 0.0); - EXPECT_EQ(dist3.pdf(0.25), 2.5); + EXPECT_DOUBLE_EQ(dist3.pdf(-0.1), 0.0); + EXPECT_DOUBLE_EQ(dist3.pdf(0.25), 2.5); } } // ppl \ No newline at end of file From 04bd73f61c689e7588a617d5485d5c9a7a2131ad Mon Sep 17 00:00:00 2001 From: Jacob Austin Date: Sat, 18 Apr 2020 14:17:39 -0400 Subject: [PATCH 18/35] added simple sampling test --- test/distribution/uniform_unittest.cpp | 11 +++++++++++ 1 file changed, 11 insertions(+) diff --git a/test/distribution/uniform_unittest.cpp b/test/distribution/uniform_unittest.cpp index fa13c601..e7366b65 100644 --- a/test/distribution/uniform_unittest.cpp +++ b/test/distribution/uniform_unittest.cpp @@ -27,4 +27,15 @@ TEST_F(uniform_dist_fixture, simple_uniform) { EXPECT_DOUBLE_EQ(dist3.pdf(0.25), 2.5); } +TEST_F(uniform_dist_fixture, uniform_sampling) { + std::random_device rd{}; + std::mt19937 gen{rd()}; + + for (int i = 0; i < 100; i++) { + double sample = dist1.sample(gen); + EXPECT_GT(sample, 0.0); + EXPECT_LT(sample, 1.0); + } +} + } // ppl \ No newline at end of file From bf549b54db05c389d4ed30aa05d1aa79bd5c6077 Mon Sep 17 00:00:00 2001 From: Jacob Austin Date: Sat, 18 Apr 2020 14:38:16 -0400 Subject: [PATCH 19/35] improved tests and distributions --- include/autoppl/distribution/normal.hpp | 23 ++++++-------------- include/autoppl/distribution/uniform.hpp | 27 ++++++++---------------- test/distribution/normal_unittest.cpp | 8 +++++++ test/distribution/uniform_unittest.cpp | 15 +++++++++++-- 4 files changed, 37 insertions(+), 36 deletions(-) diff --git a/include/autoppl/distribution/normal.hpp b/include/autoppl/distribution/normal.hpp index 8e1e6181..ed023f19 100644 --- a/include/autoppl/distribution/normal.hpp +++ b/include/autoppl/distribution/normal.hpp @@ -21,35 +21,26 @@ struct Normal { Normal(mean_type mean, var_type var) : mean_{mean}, var_{var} { - assert(static_cast(var_) > 0); + assert(this -> var() > 0); }; template value_t sample(GeneratorType& gen) const { - value_t mean, var; - mean = static_cast(mean_); - var = static_cast(var_); - - std::normal_distribution dist(mean, var); + std::normal_distribution dist(mean(), var()); return dist(gen); } dist_value_t pdf(value_t x) const { - value_t mean, var; - mean = static_cast(mean_); - var = static_cast(var_); - - return std::exp(- 0.5 * std::pow(x - mean, 2) / var) / (std::sqrt(var * 2 * M_PI)); + return std::exp(- 0.5 * std::pow(x - mean(), 2) / var()) / (std::sqrt(var() * 2 * M_PI)); } dist_value_t log_pdf(value_t x) const { - value_t mean, var; - mean = static_cast(mean_); - var = static_cast(var_); - - return (-0.5 * std::pow(x - mean, 2) / var) - 0.5 * (std::log(var) + std::log(2) + std::log(M_PI)); + return (-0.5 * std::pow(x - mean(), 2) / var()) - 0.5 * (std::log(var()) + std::log(2) + std::log(M_PI)); } + inline value_t mean() const { return static_cast(mean_);} + inline value_t var() const { return static_cast(var_);} + private: mean_type mean_; var_type var_; diff --git a/include/autoppl/distribution/uniform.hpp b/include/autoppl/distribution/uniform.hpp index 15547ace..4fb75911 100644 --- a/include/autoppl/distribution/uniform.hpp +++ b/include/autoppl/distribution/uniform.hpp @@ -16,41 +16,32 @@ struct Uniform using dist_value_t = double; Uniform(min_type min, max_type max) - : min_{min}, max_{max} { assert(static_cast(min_) < static_cast(max_)); } + : min_{min}, max_{max} { assert(this -> min() < this -> max()); } // TODO: tag this class as "TriviallySamplable"? template value_t sample(GeneratorType& gen) const { - value_t min, max; - min = static_cast(min_); - max = static_cast(max_); - - std::uniform_real_distribution dist(min, max); + std::uniform_real_distribution dist(min(), max()); return dist(gen); } dist_value_t pdf(value_t x) const { - value_t min, max; - min = static_cast(min_); - max = static_cast(max_); - - return (min < x && x < max) ? 1. / (max - min) : 0; + return (min() < x && x < max()) ? 1. / (max() - min()) : 0; } dist_value_t log_pdf(value_t x) const { - value_t min, max; - min = static_cast(min_); - max = static_cast(max_); - - return (min < x && x < max) ? - -std::log(max - min) : + return (min() < x && x < max()) ? + -std::log(max() - min()) : std::numeric_limits::lowest(); } -private: + inline value_t min() const { return static_cast(min_); } + inline value_t max() const { return static_cast(max_); } + + private: min_type min_; max_type max_; }; diff --git a/test/distribution/normal_unittest.cpp b/test/distribution/normal_unittest.cpp index 7beb42b3..939bac7d 100644 --- a/test/distribution/normal_unittest.cpp +++ b/test/distribution/normal_unittest.cpp @@ -16,6 +16,14 @@ struct normal_dist_fixture : ::testing::Test { Normal, Variable > dist2 = Normal(mu, sigma); }; +TEST_F(normal_dist_fixture, sanity_normal_test) { + EXPECT_EQ(dist1.mean(), 0.0); + EXPECT_EQ(dist1.var(), 1.0); + + EXPECT_EQ(dist2.mean(), 0.0); + EXPECT_EQ(dist2.var(), 1.0); +} + TEST_F(normal_dist_fixture, simple_gaussian) { EXPECT_DOUBLE_EQ(dist1.pdf(0.0), 0.3989422804014327); EXPECT_DOUBLE_EQ(dist1.pdf(-0.5), 0.3520653267642995); diff --git a/test/distribution/uniform_unittest.cpp b/test/distribution/uniform_unittest.cpp index e7366b65..3007a3c1 100644 --- a/test/distribution/uniform_unittest.cpp +++ b/test/distribution/uniform_unittest.cpp @@ -17,13 +17,24 @@ struct uniform_dist_fixture : ::testing::Test { Uniform, Variable > dist3 = Uniform(y, x); }; +TEST_F(uniform_dist_fixture, sanity_uniform_test) { + EXPECT_EQ(dist1.min(), 0.0); + EXPECT_EQ(dist1.max(), 1.0); + + EXPECT_EQ(dist2.min(), 0.0); + EXPECT_EQ(dist2.max(), 0.5); + + EXPECT_EQ(dist3.min(), 0.1); + EXPECT_EQ(dist3.max(), 0.5); +} + TEST_F(uniform_dist_fixture, simple_uniform) { EXPECT_DOUBLE_EQ(dist1.pdf(1.1), 0.0); + EXPECT_DOUBLE_EQ(dist1.pdf(1.0), 0.0); - EXPECT_DOUBLE_EQ(dist2.pdf(1.0), 0.0); EXPECT_DOUBLE_EQ(dist2.pdf(0.25), 2.0); + EXPECT_DOUBLE_EQ(dist2.pdf(-0.1), 0.0); - EXPECT_DOUBLE_EQ(dist3.pdf(-0.1), 0.0); EXPECT_DOUBLE_EQ(dist3.pdf(0.25), 2.5); } From 1ee259dad93e6084ce0c895c35290516c09e3029 Mon Sep 17 00:00:00 2001 From: Jacob Austin Date: Sat, 18 Apr 2020 14:56:54 -0400 Subject: [PATCH 20/35] added bernoulli distribution --- include/autoppl/distribution/bernoulli.hpp | 48 ++++++++++++++++++++++ include/autoppl/expression/model.hpp | 4 +- test/CMakeLists.txt | 1 + test/distribution/bernoulli_unittest.cpp | 42 +++++++++++++++++++ 4 files changed, 92 insertions(+), 3 deletions(-) create mode 100644 include/autoppl/distribution/bernoulli.hpp create mode 100644 test/distribution/bernoulli_unittest.cpp diff --git a/include/autoppl/distribution/bernoulli.hpp b/include/autoppl/distribution/bernoulli.hpp new file mode 100644 index 00000000..f1bc7bd1 --- /dev/null +++ b/include/autoppl/distribution/bernoulli.hpp @@ -0,0 +1,48 @@ +#pragma once +#include +#include +#include +#include + +namespace ppl { + +// TODO: change name to BernoulliDist and make class template. +// bernoulli should be a function that creates this kind of object. + +template +struct Bernoulli +{ + using value_t = int; + using dist_value_t = double; + + Bernoulli(p_type p) + : p_{p} { assert((this -> p() >= 0) && (this -> p() <= 1)); } + + template + value_t sample(GeneratorType& gen) const + { + std::bernoulli_distribution dist(p()); + return dist(gen); + } + + dist_value_t pdf(value_t x) const + { + if (x == 1) return p(); + else if (x == 0) return 1. - p(); + else return 0.0; + } + + dist_value_t log_pdf(value_t x) const + { + if (x == 1) return std::log(p()); + else if (x == 0) return std::log(1. - p()); + else return std::numeric_limits::lowest(); + } + + inline dist_value_t p() const { return static_cast(p_); } + + private: + p_type p_; +}; + +} // namespace ppl diff --git a/include/autoppl/expression/model.hpp b/include/autoppl/expression/model.hpp index a9566e8f..776f9f22 100644 --- a/include/autoppl/expression/model.hpp +++ b/include/autoppl/expression/model.hpp @@ -49,9 +49,7 @@ struct EqNode { return dist_.log_pdf(orig_var_cref_.get().get_value()); } private: - using var_cref_t = std::reference_wrapper; - using opt_var_cref_t = std::optional; - + using var_cref_t = std::reference_wrapper; var_cref_t orig_var_cref_; // (const) reference of the original var since // any configuration may be changed until right before update dist_t dist_; // distribution associated with var diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 7f93364b..670e3d9a 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -27,6 +27,7 @@ add_test(expression_unittest expression_unittest) add_executable(distribution_unittest ${CMAKE_CURRENT_SOURCE_DIR}/distribution/uniform_unittest.cpp ${CMAKE_CURRENT_SOURCE_DIR}/distribution/normal_unittest.cpp + ${CMAKE_CURRENT_SOURCE_DIR}/distribution/bernoulli_unittest.cpp ) target_compile_options(distribution_unittest PRIVATE -g -Wall -Werror -Wextra) target_include_directories(distribution_unittest PRIVATE ${GTEST_DIR}/include) diff --git a/test/distribution/bernoulli_unittest.cpp b/test/distribution/bernoulli_unittest.cpp new file mode 100644 index 00000000..d8ea2978 --- /dev/null +++ b/test/distribution/bernoulli_unittest.cpp @@ -0,0 +1,42 @@ +#include +#include + +#include +#include + +#include "gtest/gtest.h" + +namespace ppl { + +struct bernoulli_dist_fixture : ::testing::Test { +protected: + Variable x {0.6}; + + Bernoulli dist1 = Bernoulli(0.6); + Bernoulli > dist2 = Bernoulli(x); +}; + +TEST_F(bernoulli_dist_fixture, sanity_bernoulli_test) { + EXPECT_EQ(dist1.p(), 0.6); + EXPECT_EQ(dist2.p(), 0.6); +} + +TEST_F(bernoulli_dist_fixture, simple_bernoulli) { + EXPECT_DOUBLE_EQ(dist1.pdf(1), dist1.p()); + EXPECT_DOUBLE_EQ(dist1.pdf(1), 0.6); + EXPECT_DOUBLE_EQ(dist1.pdf(0), 1 - dist1.p()); + EXPECT_DOUBLE_EQ(dist1.pdf(0), 0.4); + EXPECT_DOUBLE_EQ(dist1.pdf(2), 0.0); +} + +TEST_F(bernoulli_dist_fixture, bernoulli_sampling) { + std::random_device rd{}; + std::mt19937 gen{rd()}; + + for (int i = 0; i < 100; i++) { + int sample = dist1.sample(gen); + EXPECT_TRUE(sample == 0 || sample == 1); + } +} + +} // ppl \ No newline at end of file From 67f9edfc276aab02ed930366b57c6f561448faca Mon Sep 17 00:00:00 2001 From: Jacob Austin Date: Sat, 18 Apr 2020 14:58:23 -0400 Subject: [PATCH 21/35] addressed review comments --- include/autoppl/expression/traits.hpp | 1 + test/distribution/normal_unittest.cpp | 9 +++++++++ 2 files changed, 10 insertions(+) diff --git a/include/autoppl/expression/traits.hpp b/include/autoppl/expression/traits.hpp index bb4c83d2..9d39d3d9 100644 --- a/include/autoppl/expression/traits.hpp +++ b/include/autoppl/expression/traits.hpp @@ -16,6 +16,7 @@ struct var_traits using pointer_t = typename VarType::pointer_t; using state_t = typename VarType::state_t; + // TODO may have to move this to a different class for compile-time checking static_assert(std::is_convertible_v); }; diff --git a/test/distribution/normal_unittest.cpp b/test/distribution/normal_unittest.cpp index 939bac7d..9b229a8c 100644 --- a/test/distribution/normal_unittest.cpp +++ b/test/distribution/normal_unittest.cpp @@ -32,6 +32,15 @@ TEST_F(normal_dist_fixture, simple_gaussian) { EXPECT_DOUBLE_EQ(dist1.log_pdf(0.0), std::log(dist1.pdf(0.0))); EXPECT_DOUBLE_EQ(dist1.log_pdf(-0.5), std::log(dist1.pdf(-0.5))); EXPECT_DOUBLE_EQ(dist1.log_pdf(4), std::log(dist1.pdf(4))); + + + EXPECT_DOUBLE_EQ(dist2.pdf(0.0), 0.3989422804014327); + EXPECT_DOUBLE_EQ(dist2.pdf(-0.5), 0.3520653267642995); + EXPECT_DOUBLE_EQ(dist2.pdf(4), 0.00013383022576488537); + + EXPECT_DOUBLE_EQ(dist2.log_pdf(0.0), std::log(dist2.pdf(0.0))); + EXPECT_DOUBLE_EQ(dist2.log_pdf(-0.5), std::log(dist2.pdf(-0.5))); + EXPECT_DOUBLE_EQ(dist2.log_pdf(4), std::log(dist1.pdf(4))); } } // ppl \ No newline at end of file From 43a18b8ee70d1b67c2c9697a6f185074ebbde242 Mon Sep 17 00:00:00 2001 From: James Yang Date: Sun, 19 Apr 2020 13:21:34 -0400 Subject: [PATCH 22/35] Move traits to utility directory --- include/autoppl/expression/model.hpp | 2 +- include/autoppl/{expression => utils}/traits.hpp | 0 2 files changed, 1 insertion(+), 1 deletion(-) rename include/autoppl/{expression => utils}/traits.hpp (100%) diff --git a/include/autoppl/expression/model.hpp b/include/autoppl/expression/model.hpp index 776f9f22..5f331082 100644 --- a/include/autoppl/expression/model.hpp +++ b/include/autoppl/expression/model.hpp @@ -2,7 +2,7 @@ #include #include #include -#include +#include namespace ppl { namespace details { diff --git a/include/autoppl/expression/traits.hpp b/include/autoppl/utils/traits.hpp similarity index 100% rename from include/autoppl/expression/traits.hpp rename to include/autoppl/utils/traits.hpp From 8ce9397b2fd895809fbe41bf0ca9ff2046f0a06e Mon Sep 17 00:00:00 2001 From: James Yang Date: Sun, 19 Apr 2020 13:27:57 -0400 Subject: [PATCH 23/35] Rename utils to util and remove math.h redundant dependency --- include/autoppl/distribution/normal.hpp | 2 -- include/autoppl/expression/model.hpp | 2 +- include/autoppl/{utils => util}/traits.hpp | 0 3 files changed, 1 insertion(+), 3 deletions(-) rename include/autoppl/{utils => util}/traits.hpp (100%) diff --git a/include/autoppl/distribution/normal.hpp b/include/autoppl/distribution/normal.hpp index ed023f19..14edcf14 100644 --- a/include/autoppl/distribution/normal.hpp +++ b/include/autoppl/distribution/normal.hpp @@ -1,6 +1,4 @@ #pragma once -#include - #include #include #include diff --git a/include/autoppl/expression/model.hpp b/include/autoppl/expression/model.hpp index 5f331082..baa571fa 100644 --- a/include/autoppl/expression/model.hpp +++ b/include/autoppl/expression/model.hpp @@ -2,7 +2,7 @@ #include #include #include -#include +#include namespace ppl { namespace details { diff --git a/include/autoppl/utils/traits.hpp b/include/autoppl/util/traits.hpp similarity index 100% rename from include/autoppl/utils/traits.hpp rename to include/autoppl/util/traits.hpp From 1c853ac1c997f4d1dbf01a1cf9c1be10969c1112 Mon Sep 17 00:00:00 2001 From: James Yang Date: Sun, 19 Apr 2020 13:50:06 -0400 Subject: [PATCH 24/35] Add CRTP for safety --- include/autoppl/distribution/bernoulli.hpp | 6 ++---- include/autoppl/distribution/normal.hpp | 7 +++---- include/autoppl/distribution/uniform.hpp | 6 ++---- include/autoppl/expression/model.hpp | 20 +++++++++----------- include/autoppl/expression/model_expr.hpp | 15 +++++++++++++++ include/autoppl/expression/variable.hpp | 3 ++- 6 files changed, 33 insertions(+), 24 deletions(-) create mode 100644 include/autoppl/expression/model_expr.hpp diff --git a/include/autoppl/distribution/bernoulli.hpp b/include/autoppl/distribution/bernoulli.hpp index f1bc7bd1..c4d82ad5 100644 --- a/include/autoppl/distribution/bernoulli.hpp +++ b/include/autoppl/distribution/bernoulli.hpp @@ -3,14 +3,12 @@ #include #include #include +#include namespace ppl { -// TODO: change name to BernoulliDist and make class template. -// bernoulli should be a function that creates this kind of object. - template -struct Bernoulli +struct Bernoulli : public ModelExpr> { using value_t = int; using dist_value_t = double; diff --git a/include/autoppl/distribution/normal.hpp b/include/autoppl/distribution/normal.hpp index 14edcf14..858f3927 100644 --- a/include/autoppl/distribution/normal.hpp +++ b/include/autoppl/distribution/normal.hpp @@ -3,14 +3,13 @@ #include #include #include +#include namespace ppl { -// TODO: change name to NormalDist and make class template. -// normal should be a function that creates this kind of object. - template -struct Normal { +struct Normal : public ModelExpr> +{ using value_t = double; using dist_value_t = double; diff --git a/include/autoppl/distribution/uniform.hpp b/include/autoppl/distribution/uniform.hpp index 4fb75911..16143e73 100644 --- a/include/autoppl/distribution/uniform.hpp +++ b/include/autoppl/distribution/uniform.hpp @@ -3,14 +3,12 @@ #include #include #include +#include namespace ppl { -// TODO: change name to UniformDist and make class template. -// uniform should be a function that creates this kind of object. - template -struct Uniform +struct Uniform : public ModelExpr> { using value_t = double; using dist_value_t = double; diff --git a/include/autoppl/expression/model.hpp b/include/autoppl/expression/model.hpp index baa571fa..7ba89117 100644 --- a/include/autoppl/expression/model.hpp +++ b/include/autoppl/expression/model.hpp @@ -3,6 +3,7 @@ #include #include #include +#include namespace ppl { namespace details { @@ -22,7 +23,7 @@ struct IdentityVarFunctor * that relates a var with a distribution. */ template -struct EqNode +struct EqNode : public ModelExpr> { using var_t = VarType; using dist_t = DistType; @@ -60,7 +61,7 @@ struct EqNode * "glues" two sub-model expressions. */ template -struct GlueNode +struct GlueNode : public ModelExpr> { using left_node_t = LHSNodeType; using right_node_t = RHSNodeType; @@ -98,18 +99,15 @@ struct GlueNode // Operator overloads ///////////////////////////////////////////////////////// -// TODO: all these template parameters should be constrained -// with concepts! - /* * Builds an EqNode to associate var with dist. * Ex. x |= uniform(0,1) */ template -constexpr inline auto operator|=(const VarType& var, - const DistType& dist) +constexpr inline auto operator|=(const ModelExpr& var, + const ModelExpr& dist) { - return EqNode(var, dist); + return EqNode(var.self(), dist.self()); } /* @@ -117,10 +115,10 @@ constexpr inline auto operator|=(const VarType& var, * Ex. (x |= uniform(0,1), y |= uniform(0, 2)) */ template -constexpr inline auto operator,(const LHSNodeType& lhs, - const RHSNodeType& rhs) +constexpr inline auto operator,(const ModelExpr& lhs, + const ModelExpr& rhs) { - return GlueNode(lhs, rhs); + return GlueNode(lhs.self(), rhs.self()); } } // namespace ppl diff --git a/include/autoppl/expression/model_expr.hpp b/include/autoppl/expression/model_expr.hpp new file mode 100644 index 00000000..f1f56012 --- /dev/null +++ b/include/autoppl/expression/model_expr.hpp @@ -0,0 +1,15 @@ +#pragma once + +namespace ppl { + +template +struct ModelExpr +{ + Derived& self() + { return static_cast(*this); } + + const Derived& self() const + { return static_cast(*this); } +}; + +} // namespace ppl diff --git a/include/autoppl/expression/variable.hpp b/include/autoppl/expression/variable.hpp index 3c332e9f..91021d66 100644 --- a/include/autoppl/expression/variable.hpp +++ b/include/autoppl/expression/variable.hpp @@ -1,4 +1,5 @@ #pragma once +#include namespace ppl { @@ -18,7 +19,7 @@ enum class var_state : bool { * a model expression and the users, who must supply storage of values associated with this var. */ template -struct Variable +struct Variable : public ModelExpr> { using value_t = ValueType; using pointer_t = value_t*; From fc017a51941d86250ec35c910ccb38f20484749d Mon Sep 17 00:00:00 2001 From: James Yang Date: Sun, 19 Apr 2020 13:50:29 -0400 Subject: [PATCH 25/35] Update tests for Mock classes to be ModelExpr --- test/expression/model_unittest.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/test/expression/model_unittest.cpp b/test/expression/model_unittest.cpp index 9007b3ab..24d13d2b 100644 --- a/test/expression/model_unittest.cpp +++ b/test/expression/model_unittest.cpp @@ -13,7 +13,7 @@ namespace ppl { * Mock var object for testing purposes. * Must meet some of the requirements of actual var types. */ -struct MockVar +struct MockVar : public ModelExpr { using value_t = double; using pointer_t = double*; @@ -30,7 +30,7 @@ struct MockVar * Mock distribution object for testing purposes. * Must meet some of the requirements of actual distribution types. */ -struct MockDist +struct MockDist : public ModelExpr { using value_t = double; using dist_value_t = double; From 5e92d2ef746b4694d84f47ef70b6bed3a44036ce Mon Sep 17 00:00:00 2001 From: James Yang Date: Sun, 19 Apr 2020 15:56:51 -0400 Subject: [PATCH 26/35] Work in progress (does not compile yet) but has structure --- include/autoppl/distribution/bernoulli.hpp | 4 +- include/autoppl/distribution/dist_expr.hpp | 27 ++++++ include/autoppl/distribution/normal.hpp | 4 +- include/autoppl/distribution/uniform.hpp | 4 +- include/autoppl/expr_builder.hpp | 99 ++++++++++++++++++++++ include/autoppl/expression/model.hpp | 26 ------ include/autoppl/expression/model_expr.hpp | 12 +++ include/autoppl/expression/var_expr.hpp | 27 ++++++ test/distribution/bernoulli_unittest.cpp | 2 +- test/distribution/normal_unittest.cpp | 2 +- test/distribution/uniform_unittest.cpp | 2 +- test/expression/model_unittest.cpp | 9 +- 12 files changed, 179 insertions(+), 39 deletions(-) create mode 100644 include/autoppl/distribution/dist_expr.hpp create mode 100644 include/autoppl/expr_builder.hpp create mode 100644 include/autoppl/expression/var_expr.hpp diff --git a/include/autoppl/distribution/bernoulli.hpp b/include/autoppl/distribution/bernoulli.hpp index c4d82ad5..b062988b 100644 --- a/include/autoppl/distribution/bernoulli.hpp +++ b/include/autoppl/distribution/bernoulli.hpp @@ -3,12 +3,12 @@ #include #include #include -#include +#include namespace ppl { template -struct Bernoulli : public ModelExpr> +struct Bernoulli : public DistExpr> { using value_t = int; using dist_value_t = double; diff --git a/include/autoppl/distribution/dist_expr.hpp b/include/autoppl/distribution/dist_expr.hpp new file mode 100644 index 00000000..dc559eb5 --- /dev/null +++ b/include/autoppl/distribution/dist_expr.hpp @@ -0,0 +1,27 @@ +#pragma once +#include + +namespace ppl { + +template +struct DistExpr +{ + Derived& self() + { return static_cast(*this); } + + const Derived& self() const + { return static_cast(*this); } +}; + +template +inline constexpr bool is_dist_expr_v = + std::is_convertible_v>; + +#ifdef AUTOPPL_USE_CONCEPTS +// TODO: definition should be extended with a stronger +// restriction on T with interface checking. +template +concept dist_expressable = is_dist_expr_v; +#endif + +} // namespace ppl diff --git a/include/autoppl/distribution/normal.hpp b/include/autoppl/distribution/normal.hpp index 858f3927..8e9f57db 100644 --- a/include/autoppl/distribution/normal.hpp +++ b/include/autoppl/distribution/normal.hpp @@ -3,12 +3,12 @@ #include #include #include -#include +#include namespace ppl { template -struct Normal : public ModelExpr> +struct Normal : public DistExpr> { using value_t = double; using dist_value_t = double; diff --git a/include/autoppl/distribution/uniform.hpp b/include/autoppl/distribution/uniform.hpp index 16143e73..d1a245ab 100644 --- a/include/autoppl/distribution/uniform.hpp +++ b/include/autoppl/distribution/uniform.hpp @@ -3,12 +3,12 @@ #include #include #include -#include +#include namespace ppl { template -struct Uniform : public ModelExpr> +struct Uniform : public DistExpr> { using value_t = double; using dist_value_t = double; diff --git a/include/autoppl/expr_builder.hpp b/include/autoppl/expr_builder.hpp new file mode 100644 index 00000000..e2abea32 --- /dev/null +++ b/include/autoppl/expr_builder.hpp @@ -0,0 +1,99 @@ +#pragma once +#include +#include +#include +#include + +namespace ppl { + +/* + * The purpose for these expression builders is to + * add extra type-safety and ease the user API. + */ + +//////////////////////////////////////////////////////// +// Distribution Expression Builder +//////////////////////////////////////////////////////// + +namespace details { + +using cont_param_raw_t = double; +template +inline constexpr bool is_cont_dist_param_valid = + is_var_expr_v || + std::is_convertible_v; + +} // namespace details + +#ifndef AUTOPPL_USE_CONCEPTS +/* + * Builds a Uniform expression only when the parameters + * are both valid continuous distribution parameter types. + * See var_expr.hpp for more information. + */ +template +inline constexpr auto uniform(const MinType& min_expr, + const MaxType& max_expr) +{ + static_assert(details::is_cont_dist_param_valid); + static_assert(details::is_cont_dist_param_valid); + return Uniform(min_expr, max_expr); +} +#else +#endif + +#ifndef AUTOPPL_USE_CONCEPTS +/* + * Builds a Normal expression only when the parameters + * are both valid continuous distribution parameter types. + * See var_expr.hpp for more information. + */ +template +inline constexpr auto normal(const MeanType& min_expr, + const VarianceType& max_expr) +{ + static_assert(details::is_cont_dist_param_valid); + static_assert(details::is_cont_dist_param_valid); + return Uniform(min_expr, max_expr); +} +#else +#endif + +//////////////////////////////////////////////////////// +// Model Expression Builder +//////////////////////////////////////////////////////// + +#ifndef AUTOPPL_USE_CONCEPTS +/* + * Builds an EqNode to associate var with dist + * only when var is a Variable and dist is a valid distribution expression. + * Ex. x |= uniform(0,1) + */ +template +inline constexpr auto operator|=(const Variable& var, + const DistType& dist) +{ + static_assert(is_dist_expr_v); + return EqNode(var, dist); +} +#else +#endif + +#ifndef AUTOPPL_USE_CONCEPTS +/* + * Builds a GlueNode to "glue" the left expression with the right + * only when both parameters are valid model expressions. + * Ex. (x |= uniform(0,1), y |= uniform(0, 2)) + */ +template +inline constexpr auto operator,(const LHSNodeType& lhs, + const RHSNodeType& rhs) +{ + static_assert(is_model_expr_v); + static_assert(is_model_expr_v); + return GlueNode(lhs, rhs); +} +#else +#endif + +} // namespace ppl diff --git a/include/autoppl/expression/model.hpp b/include/autoppl/expression/model.hpp index 7ba89117..4b7df09f 100644 --- a/include/autoppl/expression/model.hpp +++ b/include/autoppl/expression/model.hpp @@ -95,30 +95,4 @@ struct GlueNode : public ModelExpr> right_node_t right_node_; }; -///////////////////////////////////////////////////////// -// Operator overloads -///////////////////////////////////////////////////////// - -/* - * Builds an EqNode to associate var with dist. - * Ex. x |= uniform(0,1) - */ -template -constexpr inline auto operator|=(const ModelExpr& var, - const ModelExpr& dist) -{ - return EqNode(var.self(), dist.self()); -} - -/* - * Builds a GlueNode to "glue" the left expression with the right. - * Ex. (x |= uniform(0,1), y |= uniform(0, 2)) - */ -template -constexpr inline auto operator,(const ModelExpr& lhs, - const ModelExpr& rhs) -{ - return GlueNode(lhs.self(), rhs.self()); -} - } // namespace ppl diff --git a/include/autoppl/expression/model_expr.hpp b/include/autoppl/expression/model_expr.hpp index f1f56012..67502ae0 100644 --- a/include/autoppl/expression/model_expr.hpp +++ b/include/autoppl/expression/model_expr.hpp @@ -1,4 +1,5 @@ #pragma once +#include namespace ppl { @@ -12,4 +13,15 @@ struct ModelExpr { return static_cast(*this); } }; +template +inline constexpr bool is_model_expr_v = + std::is_convertible_v>; + +#ifdef AUTOPPL_USE_CONCEPTS +// TODO: definition should be extended with a stronger +// restriction on T with interface checking. +template +concept model_expressable = is_model_expr_v; +#endif + } // namespace ppl diff --git a/include/autoppl/expression/var_expr.hpp b/include/autoppl/expression/var_expr.hpp new file mode 100644 index 00000000..6e204243 --- /dev/null +++ b/include/autoppl/expression/var_expr.hpp @@ -0,0 +1,27 @@ +#pragma once +#include + +namespace ppl { + +template +struct VarExpr +{ + Derived& self() + { return static_cast(*this); } + + const Derived& self() const + { return static_cast(*this); } +}; + +template +inline constexpr bool is_var_expr_v = + std::is_convertible_v>; + +#ifdef AUTOPPL_USE_CONCEPTS +// TODO: definition should be extended with a stronger +// restriction on T with interface checking. +template +concept var_expressable = is_var_expr_v; +#endif + +} // namespace ppl diff --git a/test/distribution/bernoulli_unittest.cpp b/test/distribution/bernoulli_unittest.cpp index d8ea2978..47feae7a 100644 --- a/test/distribution/bernoulli_unittest.cpp +++ b/test/distribution/bernoulli_unittest.cpp @@ -39,4 +39,4 @@ TEST_F(bernoulli_dist_fixture, bernoulli_sampling) { } } -} // ppl \ No newline at end of file +} // ppl diff --git a/test/distribution/normal_unittest.cpp b/test/distribution/normal_unittest.cpp index 9b229a8c..46d4a154 100644 --- a/test/distribution/normal_unittest.cpp +++ b/test/distribution/normal_unittest.cpp @@ -43,4 +43,4 @@ TEST_F(normal_dist_fixture, simple_gaussian) { EXPECT_DOUBLE_EQ(dist2.log_pdf(4), std::log(dist1.pdf(4))); } -} // ppl \ No newline at end of file +} // ppl diff --git a/test/distribution/uniform_unittest.cpp b/test/distribution/uniform_unittest.cpp index 3007a3c1..6a426ca7 100644 --- a/test/distribution/uniform_unittest.cpp +++ b/test/distribution/uniform_unittest.cpp @@ -49,4 +49,4 @@ TEST_F(uniform_dist_fixture, uniform_sampling) { } } -} // ppl \ No newline at end of file +} // ppl diff --git a/test/expression/model_unittest.cpp b/test/expression/model_unittest.cpp index 24d13d2b..d01505b3 100644 --- a/test/expression/model_unittest.cpp +++ b/test/expression/model_unittest.cpp @@ -13,7 +13,7 @@ namespace ppl { * Mock var object for testing purposes. * Must meet some of the requirements of actual var types. */ -struct MockVar : public ModelExpr +struct MockVar { using value_t = double; using pointer_t = double*; @@ -30,7 +30,7 @@ struct MockVar : public ModelExpr * Mock distribution object for testing purposes. * Must meet some of the requirements of actual distribution types. */ -struct MockDist : public ModelExpr +struct MockDist { using value_t = double; using dist_value_t = double; @@ -51,8 +51,8 @@ struct var_dist_fixture : ::testing::Test { protected: MockVar x; - using model_t = std::decay_t; - model_t model = (x |= MockDist()); + using model_t = EqNode; + model_t model = {x, MockDist()}; double val; void reconfigure(double val) @@ -144,6 +144,7 @@ struct many_var_dist_fixture : ::testing::Test protected: MockVar x, y, z, w; double xv, yv, zv, wv; + using 2model_t = }; TEST_F(many_var_dist_fixture, two_vars) From 615d3a327d615aaf7c68b57c5fa30a2a55689931 Mon Sep 17 00:00:00 2001 From: James Yang Date: Sun, 19 Apr 2020 16:12:15 -0400 Subject: [PATCH 27/35] Fully compiles with type safety --- include/autoppl/distribution/bernoulli.hpp | 4 +-- include/autoppl/distribution/normal.hpp | 8 ++--- include/autoppl/distribution/uniform.hpp | 6 ++-- include/autoppl/expression/variable.hpp | 4 +-- test/expression/model_unittest.cpp | 34 +++++++++++++++------- 5 files changed, 34 insertions(+), 22 deletions(-) diff --git a/include/autoppl/distribution/bernoulli.hpp b/include/autoppl/distribution/bernoulli.hpp index b062988b..e034ac17 100644 --- a/include/autoppl/distribution/bernoulli.hpp +++ b/include/autoppl/distribution/bernoulli.hpp @@ -37,9 +37,9 @@ struct Bernoulli : public DistExpr> else return std::numeric_limits::lowest(); } - inline dist_value_t p() const { return static_cast(p_); } + dist_value_t p() const { return static_cast(p_); } - private: +private: p_type p_; }; diff --git a/include/autoppl/distribution/normal.hpp b/include/autoppl/distribution/normal.hpp index 8e9f57db..fc573e84 100644 --- a/include/autoppl/distribution/normal.hpp +++ b/include/autoppl/distribution/normal.hpp @@ -1,8 +1,8 @@ #pragma once #include +#include #include #include -#include #include namespace ppl { @@ -35,10 +35,10 @@ struct Normal : public DistExpr> return (-0.5 * std::pow(x - mean(), 2) / var()) - 0.5 * (std::log(var()) + std::log(2) + std::log(M_PI)); } - inline value_t mean() const { return static_cast(mean_);} - inline value_t var() const { return static_cast(var_);} + value_t mean() const { return static_cast(mean_);} + value_t var() const { return static_cast(var_);} - private: +private: mean_type mean_; var_type var_; }; diff --git a/include/autoppl/distribution/uniform.hpp b/include/autoppl/distribution/uniform.hpp index d1a245ab..1514841c 100644 --- a/include/autoppl/distribution/uniform.hpp +++ b/include/autoppl/distribution/uniform.hpp @@ -36,10 +36,10 @@ struct Uniform : public DistExpr> std::numeric_limits::lowest(); } - inline value_t min() const { return static_cast(min_); } - inline value_t max() const { return static_cast(max_); } + value_t min() const { return static_cast(min_); } + value_t max() const { return static_cast(max_); } - private: +private: min_type min_; max_type max_; }; diff --git a/include/autoppl/expression/variable.hpp b/include/autoppl/expression/variable.hpp index 91021d66..88130045 100644 --- a/include/autoppl/expression/variable.hpp +++ b/include/autoppl/expression/variable.hpp @@ -1,5 +1,5 @@ #pragma once -#include +#include namespace ppl { @@ -19,7 +19,7 @@ enum class var_state : bool { * a model expression and the users, who must supply storage of values associated with this var. */ template -struct Variable : public ModelExpr> +struct Variable : public VarExpr> { using value_t = ValueType; using pointer_t = value_t*; diff --git a/test/expression/model_unittest.cpp b/test/expression/model_unittest.cpp index d01505b3..987ceef4 100644 --- a/test/expression/model_unittest.cpp +++ b/test/expression/model_unittest.cpp @@ -144,15 +144,16 @@ struct many_var_dist_fixture : ::testing::Test protected: MockVar x, y, z, w; double xv, yv, zv, wv; - using 2model_t = + using eq_t = EqNode; }; TEST_F(many_var_dist_fixture, two_vars) { - auto model = ( - x |= MockDist(), - y |= MockDist() - ); + using model_two_t = GlueNode; + model_two_t model = { + {x, MockDist()}, + {y, MockDist()} + }; xv = 0.2; yv = 1.8; @@ -165,12 +166,23 @@ TEST_F(many_var_dist_fixture, two_vars) TEST_F(many_var_dist_fixture, four_vars) { - auto model = ( - x |= MockDist(), - y |= MockDist(), - z |= MockDist(), - w |= MockDist() - ); + using model_four_t = + GlueNode + > + >; + + model_four_t model = { + {x, MockDist()}, + { + {y, MockDist()}, + { + {z, MockDist()}, + {w, MockDist()} + } + } + }; xv = 0.2; yv = 1.8; zv = 3.2; wv = 0.3; From f09ada0eb3e5ddd6b02dd4aac304844b222f197b Mon Sep 17 00:00:00 2001 From: James Yang Date: Sun, 19 Apr 2020 16:20:20 -0400 Subject: [PATCH 28/35] Add bernoulli builder --- include/autoppl/expr_builder.hpp | 29 ++++++++++++++++++++++------- 1 file changed, 22 insertions(+), 7 deletions(-) diff --git a/include/autoppl/expr_builder.hpp b/include/autoppl/expr_builder.hpp index e2abea32..3492f1d7 100644 --- a/include/autoppl/expr_builder.hpp +++ b/include/autoppl/expr_builder.hpp @@ -17,11 +17,11 @@ namespace ppl { namespace details { -using cont_param_raw_t = double; +using cont_raw_param_t = double; template -inline constexpr bool is_cont_dist_param_valid = +inline constexpr bool is_cont_param_valid = is_var_expr_v || - std::is_convertible_v; + std::is_convertible_v; } // namespace details @@ -35,8 +35,8 @@ template inline constexpr auto uniform(const MinType& min_expr, const MaxType& max_expr) { - static_assert(details::is_cont_dist_param_valid); - static_assert(details::is_cont_dist_param_valid); + static_assert(details::is_cont_param_valid); + static_assert(details::is_cont_param_valid); return Uniform(min_expr, max_expr); } #else @@ -52,13 +52,28 @@ template inline constexpr auto normal(const MeanType& min_expr, const VarianceType& max_expr) { - static_assert(details::is_cont_dist_param_valid); - static_assert(details::is_cont_dist_param_valid); + static_assert(details::is_cont_param_valid); + static_assert(details::is_cont_param_valid); return Uniform(min_expr, max_expr); } #else #endif +#ifndef AUTOPPL_USE_CONCEPTS +/* + * Builds a Bernoulli expression only when the parameter + * is a valid discrete distribution parameter type. + * See var_expr.hpp for more information. + */ +template +inline constexpr auto bernoulli(const ProbType& p_expr) +{ + static_assert(details::is_cont_param_valid); + return Bernoulli(p_expr); +} +#else +#endif + //////////////////////////////////////////////////////// // Model Expression Builder //////////////////////////////////////////////////////// From 7341f48487cd84a331610e4d1bdf9f55ebd503f8 Mon Sep 17 00:00:00 2001 From: James Yang Date: Sun, 19 Apr 2020 16:33:30 -0400 Subject: [PATCH 29/35] Add namespace to each subdirectory --- include/autoppl/algorithm/mh.hpp | 7 +++++++ include/autoppl/distribution/bernoulli.hpp | 2 ++ include/autoppl/distribution/dist_expr.hpp | 4 +++- include/autoppl/distribution/normal.hpp | 4 +++- include/autoppl/distribution/uniform.hpp | 2 ++ include/autoppl/expr_builder.hpp | 14 +++++++++----- include/autoppl/expression/model.hpp | 2 ++ include/autoppl/expression/model_expr.hpp | 2 ++ include/autoppl/expression/var_expr.hpp | 2 ++ include/autoppl/expression/variable.hpp | 2 +- test/distribution/bernoulli_unittest.cpp | 4 +++- test/distribution/normal_unittest.cpp | 4 +++- test/distribution/uniform_unittest.cpp | 4 +++- test/expression/model_unittest.cpp | 2 ++ 14 files changed, 44 insertions(+), 11 deletions(-) create mode 100644 include/autoppl/algorithm/mh.hpp diff --git a/include/autoppl/algorithm/mh.hpp b/include/autoppl/algorithm/mh.hpp new file mode 100644 index 00000000..ffa381f9 --- /dev/null +++ b/include/autoppl/algorithm/mh.hpp @@ -0,0 +1,7 @@ +#pragma once + +namespace ppl { + + + +} // namespace ppl diff --git a/include/autoppl/distribution/bernoulli.hpp b/include/autoppl/distribution/bernoulli.hpp index e034ac17..1dfe7751 100644 --- a/include/autoppl/distribution/bernoulli.hpp +++ b/include/autoppl/distribution/bernoulli.hpp @@ -6,6 +6,7 @@ #include namespace ppl { +namespace dist { template struct Bernoulli : public DistExpr> @@ -43,4 +44,5 @@ struct Bernoulli : public DistExpr> p_type p_; }; +} // namespace dist } // namespace ppl diff --git a/include/autoppl/distribution/dist_expr.hpp b/include/autoppl/distribution/dist_expr.hpp index dc559eb5..1eb23168 100644 --- a/include/autoppl/distribution/dist_expr.hpp +++ b/include/autoppl/distribution/dist_expr.hpp @@ -2,6 +2,7 @@ #include namespace ppl { +namespace dist { template struct DistExpr @@ -15,7 +16,7 @@ struct DistExpr template inline constexpr bool is_dist_expr_v = - std::is_convertible_v>; + std::is_convertible_v>; #ifdef AUTOPPL_USE_CONCEPTS // TODO: definition should be extended with a stronger @@ -24,4 +25,5 @@ template concept dist_expressable = is_dist_expr_v; #endif +} // namespace dist } // namespace ppl diff --git a/include/autoppl/distribution/normal.hpp b/include/autoppl/distribution/normal.hpp index fc573e84..b394e803 100644 --- a/include/autoppl/distribution/normal.hpp +++ b/include/autoppl/distribution/normal.hpp @@ -6,6 +6,7 @@ #include namespace ppl { +namespace dist { template struct Normal : public DistExpr> @@ -43,4 +44,5 @@ struct Normal : public DistExpr> var_type var_; }; -} // namespace ppl +} // namespace dist +} // namespace ppl diff --git a/include/autoppl/distribution/uniform.hpp b/include/autoppl/distribution/uniform.hpp index 1514841c..4ac6f303 100644 --- a/include/autoppl/distribution/uniform.hpp +++ b/include/autoppl/distribution/uniform.hpp @@ -6,6 +6,7 @@ #include namespace ppl { +namespace dist { template struct Uniform : public DistExpr> @@ -44,4 +45,5 @@ struct Uniform : public DistExpr> max_type max_; }; +} // namespace dist } // namespace ppl diff --git a/include/autoppl/expr_builder.hpp b/include/autoppl/expr_builder.hpp index 3492f1d7..65d10ced 100644 --- a/include/autoppl/expr_builder.hpp +++ b/include/autoppl/expr_builder.hpp @@ -1,8 +1,12 @@ #pragma once -#include #include #include +#include +#include +#include #include +#include +#include namespace ppl { @@ -20,7 +24,7 @@ namespace details { using cont_raw_param_t = double; template inline constexpr bool is_cont_param_valid = - is_var_expr_v || + expr::is_var_expr_v || std::is_convertible_v; } // namespace details @@ -88,7 +92,7 @@ template inline constexpr auto operator|=(const Variable& var, const DistType& dist) { - static_assert(is_dist_expr_v); + static_assert(dist::is_dist_expr_v); return EqNode(var, dist); } #else @@ -104,8 +108,8 @@ template inline constexpr auto operator,(const LHSNodeType& lhs, const RHSNodeType& rhs) { - static_assert(is_model_expr_v); - static_assert(is_model_expr_v); + static_assert(expr::is_model_expr_v); + static_assert(expr::is_model_expr_v); return GlueNode(lhs, rhs); } #else diff --git a/include/autoppl/expression/model.hpp b/include/autoppl/expression/model.hpp index 4b7df09f..478d241a 100644 --- a/include/autoppl/expression/model.hpp +++ b/include/autoppl/expression/model.hpp @@ -6,6 +6,7 @@ #include namespace ppl { +namespace expr { namespace details { template @@ -95,4 +96,5 @@ struct GlueNode : public ModelExpr> right_node_t right_node_; }; +} // namespace expr } // namespace ppl diff --git a/include/autoppl/expression/model_expr.hpp b/include/autoppl/expression/model_expr.hpp index 67502ae0..05a26304 100644 --- a/include/autoppl/expression/model_expr.hpp +++ b/include/autoppl/expression/model_expr.hpp @@ -2,6 +2,7 @@ #include namespace ppl { +namespace expr { template struct ModelExpr @@ -24,4 +25,5 @@ template concept model_expressable = is_model_expr_v; #endif +} // namespace expr } // namespace ppl diff --git a/include/autoppl/expression/var_expr.hpp b/include/autoppl/expression/var_expr.hpp index 6e204243..c891e5b9 100644 --- a/include/autoppl/expression/var_expr.hpp +++ b/include/autoppl/expression/var_expr.hpp @@ -2,6 +2,7 @@ #include namespace ppl { +namespace expr { template struct VarExpr @@ -24,4 +25,5 @@ template concept var_expressable = is_var_expr_v; #endif +} // namespace expr } // namespace ppl diff --git a/include/autoppl/expression/variable.hpp b/include/autoppl/expression/variable.hpp index 88130045..712c8dca 100644 --- a/include/autoppl/expression/variable.hpp +++ b/include/autoppl/expression/variable.hpp @@ -19,7 +19,7 @@ enum class var_state : bool { * a model expression and the users, who must supply storage of values associated with this var. */ template -struct Variable : public VarExpr> +struct Variable : public expr::VarExpr> { using value_t = ValueType; using pointer_t = value_t*; diff --git a/test/distribution/bernoulli_unittest.cpp b/test/distribution/bernoulli_unittest.cpp index 47feae7a..267ef18e 100644 --- a/test/distribution/bernoulli_unittest.cpp +++ b/test/distribution/bernoulli_unittest.cpp @@ -7,6 +7,7 @@ #include "gtest/gtest.h" namespace ppl { +namespace dist { struct bernoulli_dist_fixture : ::testing::Test { protected: @@ -39,4 +40,5 @@ TEST_F(bernoulli_dist_fixture, bernoulli_sampling) { } } -} // ppl +} // namespace dist +} // namespace ppl diff --git a/test/distribution/normal_unittest.cpp b/test/distribution/normal_unittest.cpp index 46d4a154..36e6ee42 100644 --- a/test/distribution/normal_unittest.cpp +++ b/test/distribution/normal_unittest.cpp @@ -7,6 +7,7 @@ #include "gtest/gtest.h" namespace ppl { +namespace dist { struct normal_dist_fixture : ::testing::Test { protected: @@ -43,4 +44,5 @@ TEST_F(normal_dist_fixture, simple_gaussian) { EXPECT_DOUBLE_EQ(dist2.log_pdf(4), std::log(dist1.pdf(4))); } -} // ppl +} // namespace dist +} // namespace ppl diff --git a/test/distribution/uniform_unittest.cpp b/test/distribution/uniform_unittest.cpp index 6a426ca7..cf1d7b4b 100644 --- a/test/distribution/uniform_unittest.cpp +++ b/test/distribution/uniform_unittest.cpp @@ -7,6 +7,7 @@ #include "gtest/gtest.h" namespace ppl { +namespace dist { struct uniform_dist_fixture : ::testing::Test { protected: @@ -49,4 +50,5 @@ TEST_F(uniform_dist_fixture, uniform_sampling) { } } -} // ppl +} // namespace dist +} // namespace ppl diff --git a/test/expression/model_unittest.cpp b/test/expression/model_unittest.cpp index 987ceef4..df3a6f6a 100644 --- a/test/expression/model_unittest.cpp +++ b/test/expression/model_unittest.cpp @@ -4,6 +4,7 @@ #include namespace ppl { +namespace expr { ////////////////////////////////////////////////////// // Model with one RV TESTS @@ -196,4 +197,5 @@ TEST_F(many_var_dist_fixture, four_vars) + std::log(zv) + std::log(wv)); } +} // namespace expr } // namespace ppl From 2db24c78e6c0508d96eda8eb251fd77dcfed956e Mon Sep 17 00:00:00 2001 From: James Yang Date: Sun, 19 Apr 2020 22:44:54 -0400 Subject: [PATCH 30/35] Add metropolis hastings! --- include/autoppl/algorithm/mh.hpp | 131 +++++++++++++++++++++ include/autoppl/distribution/bernoulli.hpp | 21 ++-- include/autoppl/distribution/density.hpp | 90 ++++++++++++++ include/autoppl/distribution/normal.hpp | 41 +++---- include/autoppl/distribution/uniform.hpp | 24 ++-- include/autoppl/expr_builder.hpp | 20 ++-- include/autoppl/expression/model.hpp | 42 +++++-- include/autoppl/util/traits.hpp | 12 +- test/CMakeLists.txt | 15 +++ test/algorithm/mh_unittest.cpp | 74 ++++++++++++ test/distribution/normal_unittest.cpp | 4 +- test/expression/model_unittest.cpp | 78 ++++++++---- 12 files changed, 465 insertions(+), 87 deletions(-) create mode 100644 include/autoppl/distribution/density.hpp create mode 100644 test/algorithm/mh_unittest.cpp diff --git a/include/autoppl/algorithm/mh.hpp b/include/autoppl/algorithm/mh.hpp index ffa381f9..5f59b961 100644 --- a/include/autoppl/algorithm/mh.hpp +++ b/include/autoppl/algorithm/mh.hpp @@ -1,7 +1,138 @@ #pragma once +#include +#include +#include +#include +#include +#include + +/* + * Assumptions: + * - every variable referenced in model is of type Variable + */ namespace ppl { +namespace details { + +struct MHData +{ + double next; + // TODO: maybe keep an array for batch sampling? +}; + +} // namespace details + +/* + * Metropolis-Hastings algorithm to sample from posterior distribution. + * The posterior distribution is a constant multiple of model.pdf(). + * Any variables that model references which are in state "parameter" + * is sampled and in state "data" are not. + * So, model.pdf() is proportional to p(parameters... | data...). + * + * User must ensure that they allocated at least as large as n_sample + * in the storage associated with every parameter referenced in model. + */ +template +inline void mh_posterior(ModelType& model, + double n_sample, + double stddev = 1.0, + double seed = 0 + ) +{ + using data_t = details::MHData; + + // set-up auxiliary tools + std::mt19937 gen(seed); + std::uniform_real_distribution unif_sampler(0., 1.); + + // get number of parameters to sample + size_t n_params = 0.; + auto get_n_params = [&](auto& eq_node) { + auto& var = eq_node.get_variable(); + using var_t = std::decay_t; + using state_t = typename var_traits::state_t; + n_params += (var.get_state() == state_t::parameter); + }; + model.traverse(get_n_params); + + // vector of parameter-related data with candidate + std::vector params(n_params); + double curr_log_pdf = model.log_pdf(); + auto params_it = params.begin(); + + for (size_t iter = 0; iter < n_sample; ++iter) { + + double log_alpha = -curr_log_pdf; + + // generate next candidates and place them in parameter + // variables as next values; update log_alpha + // The old values are temporary stored in the params vector. + auto get_candidate = [=, &gen](auto& eq_node) mutable { + auto& var = eq_node.get_variable(); + using var_t = std::decay_t; + using value_t = typename var_traits::value_t; + using state_t = typename var_traits::state_t; + + if (var.get_state() == state_t::parameter) { + auto curr = static_cast(var); + std::normal_distribution norm_sampler(curr, stddev); + + // sample new candidate, place old value in params, + // fill next candidate in var, and update log_alpha + auto cand = norm_sampler(gen); + params_it->next = curr; + var.set_value(cand); + log_alpha -= dist::NormalBase::log_pdf(cand, curr, stddev); + log_alpha += dist::NormalBase::log_pdf(curr, cand, stddev); + + ++params_it; + } + }; + model.traverse(get_candidate); + + // compute next candidate log pdf and update log_alpha + double cand_log_pdf = model.log_pdf(); + log_alpha += cand_log_pdf; + + // accept + if (std::log(unif_sampler(gen)) <= log_alpha) { + // "current" sample for next iteration is already in the variables + // simply append to storage + auto add_to_storage = [iter](auto& eq_node) { + auto& var = eq_node.get_variable(); + using var_t = std::decay_t; + using state_t = typename var_traits::state_t; + if (var.get_state() == state_t::parameter) { + auto storage = var.get_storage(); + storage[iter] = var.get_value(); + } + }; + model.traverse(add_to_storage); + + // update current log pdf for next iteration + curr_log_pdf = cand_log_pdf; + } + // reject + else { + // "current" sample for next iteration must be moved back from + // params vector into variables. + // Then, simply append to storage + auto add_to_storage = [params_it, iter](auto& eq_node) mutable { + auto& var = eq_node.get_variable(); + using var_t = std::decay_t; + using state_t = typename var_traits::state_t; + if (var.get_state() == state_t::parameter) { + var.set_value(params_it->next); + auto storage = var.get_storage(); + storage[iter] = var.get_value(); + ++params_it; + } + }; + model.traverse(add_to_storage); + } + } +} } // namespace ppl diff --git a/include/autoppl/distribution/bernoulli.hpp b/include/autoppl/distribution/bernoulli.hpp index 1dfe7751..18424b79 100644 --- a/include/autoppl/distribution/bernoulli.hpp +++ b/include/autoppl/distribution/bernoulli.hpp @@ -3,7 +3,9 @@ #include #include #include +#include #include +#include namespace ppl { namespace dist { @@ -11,8 +13,9 @@ namespace dist { template struct Bernoulli : public DistExpr> { - using value_t = int; - using dist_value_t = double; + using value_t = uint8_t; + using param_value_t = typename var_traits::value_t; + using dist_value_t = typename BernoulliBase::dist_value_t; Bernoulli(p_type p) : p_{p} { assert((this -> p() >= 0) && (this -> p() <= 1)); } @@ -25,20 +28,12 @@ struct Bernoulli : public DistExpr> } dist_value_t pdf(value_t x) const - { - if (x == 1) return p(); - else if (x == 0) return 1. - p(); - else return 0.0; - } + { return BernoulliBase::pdf(x, p()); } dist_value_t log_pdf(value_t x) const - { - if (x == 1) return std::log(p()); - else if (x == 0) return std::log(1. - p()); - else return std::numeric_limits::lowest(); - } + { return BernoulliBase::log_pdf(x, p()); } - dist_value_t p() const { return static_cast(p_); } + param_value_t p() const { return static_cast(p_); } private: p_type p_; diff --git a/include/autoppl/distribution/density.hpp b/include/autoppl/distribution/density.hpp new file mode 100644 index 00000000..492d0827 --- /dev/null +++ b/include/autoppl/distribution/density.hpp @@ -0,0 +1,90 @@ +#pragma once +#include + +namespace ppl { +namespace dist { + +/* + * The Base objects contain static member functions + * that compute pdf and log_pdf. + * These are useful stand-alone functions and the distribution objects + * such as Uniform and Normal simply wrap these functions. + */ + +/* + * Continuous Distributions + */ + +struct UniformBase +{ + using dist_value_t = double; + + template + static dist_value_t pdf(ValueType x, + ParamValueType min, + ParamValueType max) + { + return (min < x && x < max) ? 1. / (max - min) : 0; + } + + template + static dist_value_t log_pdf(ValueType x, + ParamValueType min, + ParamValueType max) + { + return (min < x && x < max) ? + -std::log(max - min) : + std::numeric_limits::lowest(); + } +}; + +struct NormalBase +{ + using dist_value_t = double; + + template + static dist_value_t pdf(ValueType x, + ParamValueType mean, + ParamValueType stddev) + { + dist_value_t z_score = (x - mean) / stddev; + return std::exp(- 0.5 * z_score * z_score) / (stddev * std::sqrt(2 * M_PI)); + } + + template + static dist_value_t log_pdf(ValueType x, + ParamValueType mean, + ParamValueType stddev) + { + dist_value_t z_score = (x - mean) / stddev; + return -0.5 * ((z_score * z_score) + std::log(stddev * stddev * 2 * M_PI)); + } +}; + +/* + * Discrete Distributions + */ + +struct BernoulliBase +{ + using dist_value_t = double; + + template + static dist_value_t pdf(ValueType x, ParamValueType p) + { + if (x == 1) return p; + else if (x == 0) return 1. - p; + else return 0.0; + } + + template + static dist_value_t log_pdf(ValueType x, ParamValueType p) + { + if (x == 1) return std::log(p); + else if (x == 0) return std::log(1. - p); + else return std::numeric_limits::lowest(); + } +}; + +} // namespace dist +} // namespace ppl diff --git a/include/autoppl/distribution/normal.hpp b/include/autoppl/distribution/normal.hpp index b394e803..f069e693 100644 --- a/include/autoppl/distribution/normal.hpp +++ b/include/autoppl/distribution/normal.hpp @@ -3,45 +3,46 @@ #include #include #include +#include #include +#include namespace ppl { namespace dist { -template -struct Normal : public DistExpr> +template +struct Normal : public DistExpr> { using value_t = double; - using dist_value_t = double; - - static_assert(std::is_convertible_v); - static_assert(std::is_convertible_v); - - Normal(mean_type mean, var_type var) - : mean_{mean}, var_{var} { - assert(this -> var() > 0); + using param_value_t = std::common_type_t< + typename var_traits::value_t, + typename var_traits::value_t + >; + using dist_value_t = typename NormalBase::dist_value_t; + + Normal(mean_type mean, stddev_type stddev) + : mean_{mean}, stddev_{stddev} { + assert(this -> stddev() > 0); }; template value_t sample(GeneratorType& gen) const { - std::normal_distribution dist(mean(), var()); + std::normal_distribution dist(mean(), stddev()); return dist(gen); } - dist_value_t pdf(value_t x) const { - return std::exp(- 0.5 * std::pow(x - mean(), 2) / var()) / (std::sqrt(var() * 2 * M_PI)); - } + dist_value_t pdf(value_t x) const + { return NormalBase::pdf(x, mean(), stddev()); } - dist_value_t log_pdf(value_t x) const { - return (-0.5 * std::pow(x - mean(), 2) / var()) - 0.5 * (std::log(var()) + std::log(2) + std::log(M_PI)); - } + dist_value_t log_pdf(value_t x) const + { return NormalBase::log_pdf(x, mean(), stddev()); } - value_t mean() const { return static_cast(mean_);} - value_t var() const { return static_cast(var_);} + param_value_t mean() const { return static_cast(mean_);} + param_value_t stddev() const { return static_cast(stddev_);} private: mean_type mean_; - var_type var_; + stddev_type stddev_; }; } // namespace dist diff --git a/include/autoppl/distribution/uniform.hpp b/include/autoppl/distribution/uniform.hpp index 4ac6f303..41c830c3 100644 --- a/include/autoppl/distribution/uniform.hpp +++ b/include/autoppl/distribution/uniform.hpp @@ -3,7 +3,9 @@ #include #include #include +#include #include +#include namespace ppl { namespace dist { @@ -12,7 +14,11 @@ template struct Uniform : public DistExpr> { using value_t = double; - using dist_value_t = double; + using param_value_t = std::common_type_t< + typename var_traits::value_t, + typename var_traits::value_t + >; + using dist_value_t = typename UniformBase::dist_value_t; Uniform(min_type min, max_type max) : min_{min}, max_{max} { assert(this -> min() < this -> max()); } @@ -21,24 +27,18 @@ struct Uniform : public DistExpr> template value_t sample(GeneratorType& gen) const { - std::uniform_real_distribution dist(min(), max()); + std::uniform_real_distribution dist(min(), max()); return dist(gen); } dist_value_t pdf(value_t x) const - { - return (min() < x && x < max()) ? 1. / (max() - min()) : 0; - } + { return UniformBase::pdf(x, min(), max()); } dist_value_t log_pdf(value_t x) const - { - return (min() < x && x < max()) ? - -std::log(max() - min()) : - std::numeric_limits::lowest(); - } + { return UniformBase::log_pdf(x, min(), max()); } - value_t min() const { return static_cast(min_); } - value_t max() const { return static_cast(max_); } + param_value_t min() const { return static_cast(min_); } + param_value_t max() const { return static_cast(max_); } private: min_type min_; diff --git a/include/autoppl/expr_builder.hpp b/include/autoppl/expr_builder.hpp index 65d10ced..8fd0964f 100644 --- a/include/autoppl/expr_builder.hpp +++ b/include/autoppl/expr_builder.hpp @@ -41,7 +41,7 @@ inline constexpr auto uniform(const MinType& min_expr, { static_assert(details::is_cont_param_valid); static_assert(details::is_cont_param_valid); - return Uniform(min_expr, max_expr); + return dist::Uniform(min_expr, max_expr); } #else #endif @@ -52,13 +52,13 @@ inline constexpr auto uniform(const MinType& min_expr, * are both valid continuous distribution parameter types. * See var_expr.hpp for more information. */ -template -inline constexpr auto normal(const MeanType& min_expr, - const VarianceType& max_expr) +template +inline constexpr auto normal(const MeanType& mean_expr, + const StddevType& stddev_expr) { static_assert(details::is_cont_param_valid); - static_assert(details::is_cont_param_valid); - return Uniform(min_expr, max_expr); + static_assert(details::is_cont_param_valid); + return dist::Normal(mean_expr, stddev_expr); } #else #endif @@ -73,7 +73,7 @@ template inline constexpr auto bernoulli(const ProbType& p_expr) { static_assert(details::is_cont_param_valid); - return Bernoulli(p_expr); + return dist::Bernoulli(p_expr); } #else #endif @@ -89,11 +89,11 @@ inline constexpr auto bernoulli(const ProbType& p_expr) * Ex. x |= uniform(0,1) */ template -inline constexpr auto operator|=(const Variable& var, +inline constexpr auto operator|=(Variable& var, const DistType& dist) { static_assert(dist::is_dist_expr_v); - return EqNode(var, dist); + return expr::EqNode(var, dist); } #else #endif @@ -110,7 +110,7 @@ inline constexpr auto operator,(const LHSNodeType& lhs, { static_assert(expr::is_model_expr_v); static_assert(expr::is_model_expr_v); - return GlueNode(lhs, rhs); + return expr::GlueNode(lhs, rhs); } #else #endif diff --git a/include/autoppl/expression/model.hpp b/include/autoppl/expression/model.hpp index 478d241a..5e06c6ad 100644 --- a/include/autoppl/expression/model.hpp +++ b/include/autoppl/expression/model.hpp @@ -30,31 +30,46 @@ struct EqNode : public ModelExpr> using dist_t = DistType; using dist_value_t = typename dist_traits::dist_value_t; - EqNode(const var_t& var, + EqNode(var_t& var, const dist_t& dist) noexcept - : orig_var_cref_{var} + : orig_var_ref_{var} , dist_{dist} {} + /* + * Generic traversal function. + * Assumes that eq_f is a function that takes in 1 parameter, + * which is simply the current object. + */ + template + void traverse(EqNodeFunc&& eq_f) + { + using this_t = EqNode; + eq_f(static_cast(*this)); + } + /* * Compute pdf of underlying distribution with underlying value. * Assumes that underlying value has been assigned properly. */ dist_value_t pdf() const - { return dist_.pdf(orig_var_cref_.get().get_value()); } + { return dist_.pdf(orig_var_ref_.get().get_value()); } /* * Compute log-pdf of underlying distribution with underlying value. * Assumes that underlying value has been assigned properly. */ dist_value_t log_pdf() const - { return dist_.log_pdf(orig_var_cref_.get().get_value()); } + { return dist_.log_pdf(orig_var_ref_.get().get_value()); } + + auto& get_variable() { return orig_var_ref_.get(); } + const auto& get_distribution() const { return dist_; } private: - using var_cref_t = std::reference_wrapper; - var_cref_t orig_var_cref_; // (const) reference of the original var since - // any configuration may be changed until right before update - dist_t dist_; // distribution associated with var + using var_ref_t = std::reference_wrapper; + var_ref_t orig_var_ref_; // reference of the original var since + // any configuration may be changed until right before update + dist_t dist_; // distribution associated with var }; /* @@ -77,6 +92,17 @@ struct GlueNode : public ModelExpr> , right_node_{rhs} {} + /* + * Generic traversal function. + * Recursively traverses left then right. + */ + template + void traverse(EqNodeFunc&& eq_f) + { + left_node_.traverse(eq_f); + right_node_.traverse(eq_f); + } + /* * Computes left node joint pdf then right node joint pdf * and returns the product of the two. diff --git a/include/autoppl/util/traits.hpp b/include/autoppl/util/traits.hpp index 9d39d3d9..5e812978 100644 --- a/include/autoppl/util/traits.hpp +++ b/include/autoppl/util/traits.hpp @@ -20,11 +20,19 @@ struct var_traits static_assert(std::is_convertible_v); }; +// Specialization: when double or int, considered "trivial" variable. +// TODO: this was a quick fix for generalizing distribution value_t. +template <> +struct var_traits +{ + using value_t = double; +}; + template struct dist_traits { - using value_t = typename DistType::value_t; - using dist_value_t = typename DistType::dist_value_t; + using value_t = typename DistType::value_t; // generated value type + using dist_value_t = typename DistType::dist_value_t; // pdf value type }; template diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 670e3d9a..0383af60 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -36,3 +36,18 @@ if (AUTOPPL_ENABLE_TEST_COVERAGE) endif() target_link_libraries(distribution_unittest gtest_main pthread ${PROJECT_NAME}) add_test(distribution_unittest distribution_unittest) + +###################################################### +# Algorithm Test +###################################################### + +add_executable(algorithm_unittest + ${CMAKE_CURRENT_SOURCE_DIR}/algorithm/mh_unittest.cpp + ) +target_compile_options(algorithm_unittest PRIVATE -g -Wall -Werror -Wextra) +target_include_directories(algorithm_unittest PRIVATE ${GTEST_DIR}/include) +if (AUTOPPL_ENABLE_TEST_COVERAGE) + target_link_libraries(algorithm_unittest gcov) +endif() +target_link_libraries(algorithm_unittest gtest_main pthread ${PROJECT_NAME}) +add_test(algorithm_unittest algorithm_unittest) diff --git a/test/algorithm/mh_unittest.cpp b/test/algorithm/mh_unittest.cpp new file mode 100644 index 00000000..c1299d10 --- /dev/null +++ b/test/algorithm/mh_unittest.cpp @@ -0,0 +1,74 @@ +#include "gtest/gtest.h" +#include +#include +#include +#include + +namespace ppl { + +template +inline void plot_hist(const ArrayType& arr, + double step_size = .5) +{ + constexpr size_t nstars = 100; // maximum number of stars to distribute + + const int64_t min = std::floor(*std::min_element(arr.begin(), arr.end())); + const int64_t max = std::ceil(*std::max_element(arr.begin(), arr.end())); + const uint64_t range = max - min; + const uint64_t n_hist = std::floor(range/step_size); + + // keeps count for each histogram bar + std::vector counter(n_hist, 0); + + for (auto x : arr) { + ++counter[std::floor((x - min) / step_size)]; + } + + EXPECT_EQ(std::accumulate(counter.begin(), counter.end(), 0), (int) arr.size()); + + for (size_t i = 0; i < n_hist; ++i) { + std::cout << i << "-" << (i+1) << ": " << '\t'; + std::cout << std::string(counter[i] * nstars/arr.size(), '*') << std::endl; + } + +} + +TEST(mh_unittest, plot_hist_sanity) +{ + static constexpr size_t sample_size = 5000; + std::array storage = {0.}; + std::normal_distribution normal_sampler(0., 1.); + std::mt19937 gen; + for (size_t i = 0; i < sample_size; ++i) { + storage[i] = normal_sampler(gen); + } + plot_hist(storage); +} + +struct mh_fixture : ::testing::Test +{ +protected: + static constexpr size_t sample_size = 5000; + std::array storage = {0.}; + Variable theta; + + mh_fixture() + : theta{storage.data()} + {} +}; + +TEST_F(mh_fixture, sample_std_normal) +{ + auto model = (theta |= normal(0., 1.)); + mh_posterior(model, sample_size); + plot_hist(storage); +} + +TEST_F(mh_fixture, sample_uniform) +{ + auto model = (theta |= uniform(0., 1.)); + mh_posterior(model, sample_size); + plot_hist(storage, 0.05); +} + +} // namespace ppl diff --git a/test/distribution/normal_unittest.cpp b/test/distribution/normal_unittest.cpp index 36e6ee42..3dea0eff 100644 --- a/test/distribution/normal_unittest.cpp +++ b/test/distribution/normal_unittest.cpp @@ -19,10 +19,10 @@ struct normal_dist_fixture : ::testing::Test { TEST_F(normal_dist_fixture, sanity_normal_test) { EXPECT_EQ(dist1.mean(), 0.0); - EXPECT_EQ(dist1.var(), 1.0); + EXPECT_EQ(dist1.stddev(), 1.0); EXPECT_EQ(dist2.mean(), 0.0); - EXPECT_EQ(dist2.var(), 1.0); + EXPECT_EQ(dist2.stddev(), 1.0); } TEST_F(normal_dist_fixture, simple_gaussian) { diff --git a/test/expression/model_unittest.cpp b/test/expression/model_unittest.cpp index df3a6f6a..4cc39eec 100644 --- a/test/expression/model_unittest.cpp +++ b/test/expression/model_unittest.cpp @@ -10,6 +10,14 @@ namespace expr { // Model with one RV TESTS ////////////////////////////////////////////////////// +/* + * Mock state class for testing purposes. + */ +enum class MockState { + data, + parameter +}; + /* * Mock var object for testing purposes. * Must meet some of the requirements of actual var types. @@ -18,13 +26,19 @@ struct MockVar { using value_t = double; using pointer_t = double*; - using state_t = void; + using state_t = MockState; void set_value(double val) { value_ = val; } - double get_value() const { return value_; } + value_t get_value() const { return value_; } + + void set_state(state_t state) { state_ = state; } + state_t get_state() const { return state_; } + + operator value_t() { return value_; } private: double value_; + state_t state_ = state_t::parameter; }; /* @@ -146,27 +160,14 @@ struct many_var_dist_fixture : ::testing::Test MockVar x, y, z, w; double xv, yv, zv, wv; using eq_t = EqNode; -}; -TEST_F(many_var_dist_fixture, two_vars) -{ using model_two_t = GlueNode; - model_two_t model = { + model_two_t model_two = { {x, MockDist()}, {y, MockDist()} - }; - xv = 0.2; yv = 1.8; - - x.set_value(xv); - y.set_value(yv); - - EXPECT_EQ(model.pdf(), xv * yv); - EXPECT_EQ(model.log_pdf(), std::log(xv) + std::log(yv)); -} + }; -TEST_F(many_var_dist_fixture, four_vars) -{ using model_four_t = GlueNode >; - model_four_t model = { + model_four_t model_four = { {x, MockDist()}, { {y, MockDist()}, @@ -184,7 +185,21 @@ TEST_F(many_var_dist_fixture, four_vars) } } }; +}; + +TEST_F(many_var_dist_fixture, two_vars_pdf) +{ + xv = 0.2; yv = 1.8; + + x.set_value(xv); + y.set_value(yv); + + EXPECT_EQ(model_two.pdf(), xv * yv); + EXPECT_EQ(model_two.log_pdf(), std::log(xv) + std::log(yv)); +} +TEST_F(many_var_dist_fixture, four_vars_pdf) +{ xv = 0.2; yv = 1.8; zv = 3.2; wv = 0.3; x.set_value(xv); @@ -192,10 +207,33 @@ TEST_F(many_var_dist_fixture, four_vars) z.set_value(zv); w.set_value(wv); - EXPECT_EQ(model.pdf(), xv * yv * zv * wv); - EXPECT_EQ(model.log_pdf(), std::log(xv) + std::log(yv) + EXPECT_EQ(model_four.pdf(), xv * yv * zv * wv); + EXPECT_EQ(model_four.log_pdf(), std::log(xv) + std::log(yv) + std::log(zv) + std::log(wv)); } +TEST_F(many_var_dist_fixture, four_vars_traverse_count_params) +{ + int count = 0; + z.set_state(MockState::data); + model_four.traverse([&](auto& model) { + using var_t = std::decay_t; + using state_t = typename var_traits::state_t; + count += (model.get_variable().get_state() == state_t::parameter); + }); + EXPECT_EQ(count, 3); +} + +TEST_F(many_var_dist_fixture, four_vars_traverse_pdf) +{ + double actual = 1.; + model_four.traverse([&](auto& model) { + auto& var = model.get_variable(); + auto& dist = model.get_distribution(); + actual *= dist.pdf(var); + }); + EXPECT_EQ(actual, model_four.pdf()); +} + } // namespace expr } // namespace ppl From 548f5de7b6b2804c79570f6d49ccd512d0417727 Mon Sep 17 00:00:00 2001 From: James Yang Date: Sun, 19 Apr 2020 23:06:29 -0400 Subject: [PATCH 31/35] Set value at constructor for variable sets to data --- include/autoppl/expression/variable.hpp | 11 ++++++----- 1 file changed, 6 insertions(+), 5 deletions(-) diff --git a/include/autoppl/expression/variable.hpp b/include/autoppl/expression/variable.hpp index 712c8dca..eeca0092 100644 --- a/include/autoppl/expression/variable.hpp +++ b/include/autoppl/expression/variable.hpp @@ -28,21 +28,22 @@ struct Variable : public expr::VarExpr> // constructors Variable(value_t value, - pointer_t storage_ptr) noexcept + pointer_t storage_ptr, + state_t state) noexcept : value_{value} , storage_ptr_{storage_ptr} - , state_{state_t::parameter} + , state_{state} {} Variable(pointer_t storage_ptr) noexcept - : Variable(0, storage_ptr) + : Variable(0, storage_ptr, state_t::parameter) {} Variable(value_t value) noexcept - : Variable(value, nullptr) {} + : Variable(value, nullptr, state_t::data) {} Variable() noexcept - : Variable(0, nullptr) + : Variable(0, nullptr, state_t::parameter) {} void set_value(value_t value) { value_ = value; } From 20c486fff62ef6d4289848bf5f5eedec76551fd2 Mon Sep 17 00:00:00 2001 From: James Yang Date: Mon, 20 Apr 2020 17:45:58 -0400 Subject: [PATCH 32/35] Fully refactored and conceptualized --- include/autoppl/algorithm/mh.hpp | 67 +++----- include/autoppl/distribution/dist_expr.hpp | 29 ---- include/autoppl/expr_builder.hpp | 114 +++++++++---- .../distribution/bernoulli.hpp | 20 +-- .../{ => expression}/distribution/density.hpp | 5 +- .../{ => expression}/distribution/normal.hpp | 23 +-- .../{ => expression}/distribution/uniform.hpp | 23 +-- .../autoppl/expression/{ => model}/model.hpp | 32 ++-- include/autoppl/expression/model_expr.hpp | 29 ---- include/autoppl/expression/var_expr.hpp | 29 ---- .../autoppl/expression/variable/constant.hpp | 20 +++ .../expression/variable/variable_viewer.hpp | 33 ++++ include/autoppl/util/concept.hpp | 127 ++++++++++++++ include/autoppl/util/dist_expr_traits.hpp | 60 +++++++ include/autoppl/util/model_expr_traits.hpp | 43 +++++ include/autoppl/util/traits.hpp | 40 +---- include/autoppl/util/var_expr_traits.hpp | 47 ++++++ include/autoppl/util/var_traits.hpp | 43 +++++ include/autoppl/{expression => }/variable.hpp | 9 +- test/CMakeLists.txt | 106 +++++++++--- test/algorithm/mh_unittest.cpp | 64 +++---- test/distribution/bernoulli_unittest.cpp | 44 ----- test/distribution/normal_unittest.cpp | 48 ------ test/distribution/uniform_unittest.cpp | 54 ------ test/expr_builder_unittest.cpp | 55 ++++++ .../distribution/bernoulli_unittest.cpp | 65 ++++++++ .../distribution/density_unittest.cpp | 157 ++++++++++++++++++ .../distribution/normal_unittest.cpp | 68 ++++++++ .../distribution/uniform_unittest.cpp | 70 ++++++++ .../expression/{ => model}/model_unittest.cpp | 144 ++++------------ .../expression/variable/constant_unittest.cpp | 30 ++++ .../variable/variable_viewer_unittest.cpp | 32 ++++ test/testutil/mock_types.hpp | 118 +++++++++++++ test/testutil/sample_tools.hpp | 53 ++++++ test/util/concept_unittest.cpp | 91 ++++++++++ test/util/dist_expr_traits_unittest.cpp | 25 +++ test/util/var_expr_traits_unittest.cpp | 24 +++ test/util/var_traits_unittest.cpp | 24 +++ 38 files changed, 1500 insertions(+), 565 deletions(-) delete mode 100644 include/autoppl/distribution/dist_expr.hpp rename include/autoppl/{ => expression}/distribution/bernoulli.hpp (66%) rename include/autoppl/{ => expression}/distribution/density.hpp (97%) rename include/autoppl/{ => expression}/distribution/normal.hpp (67%) rename include/autoppl/{ => expression}/distribution/uniform.hpp (68%) rename include/autoppl/expression/{ => model}/model.hpp (82%) delete mode 100644 include/autoppl/expression/model_expr.hpp delete mode 100644 include/autoppl/expression/var_expr.hpp create mode 100644 include/autoppl/expression/variable/constant.hpp create mode 100644 include/autoppl/expression/variable/variable_viewer.hpp create mode 100644 include/autoppl/util/concept.hpp create mode 100644 include/autoppl/util/dist_expr_traits.hpp create mode 100644 include/autoppl/util/model_expr_traits.hpp create mode 100644 include/autoppl/util/var_expr_traits.hpp create mode 100644 include/autoppl/util/var_traits.hpp rename include/autoppl/{expression => }/variable.hpp (90%) delete mode 100644 test/distribution/bernoulli_unittest.cpp delete mode 100644 test/distribution/normal_unittest.cpp delete mode 100644 test/distribution/uniform_unittest.cpp create mode 100644 test/expr_builder_unittest.cpp create mode 100644 test/expression/distribution/bernoulli_unittest.cpp create mode 100644 test/expression/distribution/density_unittest.cpp create mode 100644 test/expression/distribution/normal_unittest.cpp create mode 100644 test/expression/distribution/uniform_unittest.cpp rename test/expression/{ => model}/model_unittest.cpp (53%) create mode 100644 test/expression/variable/constant_unittest.cpp create mode 100644 test/expression/variable/variable_viewer_unittest.cpp create mode 100644 test/testutil/mock_types.hpp create mode 100644 test/testutil/sample_tools.hpp create mode 100644 test/util/concept_unittest.cpp create mode 100644 test/util/dist_expr_traits_unittest.cpp create mode 100644 test/util/var_expr_traits_unittest.cpp create mode 100644 test/util/var_traits_unittest.cpp diff --git a/include/autoppl/algorithm/mh.hpp b/include/autoppl/algorithm/mh.hpp index 5f59b961..1a83d3f8 100644 --- a/include/autoppl/algorithm/mh.hpp +++ b/include/autoppl/algorithm/mh.hpp @@ -3,8 +3,8 @@ #include #include #include -#include -#include +#include +#include /* * Assumptions: @@ -50,7 +50,7 @@ inline void mh_posterior(ModelType& model, auto get_n_params = [&](auto& eq_node) { auto& var = eq_node.get_variable(); using var_t = std::decay_t; - using state_t = typename var_traits::state_t; + using state_t = typename util::var_traits::state_t; n_params += (var.get_state() == state_t::parameter); }; model.traverse(get_n_params); @@ -70,11 +70,10 @@ inline void mh_posterior(ModelType& model, auto get_candidate = [=, &gen](auto& eq_node) mutable { auto& var = eq_node.get_variable(); using var_t = std::decay_t; - using value_t = typename var_traits::value_t; - using state_t = typename var_traits::state_t; + using state_t = typename util::var_traits::state_t; if (var.get_state() == state_t::parameter) { - auto curr = static_cast(var); + auto curr = var.get_value(); std::normal_distribution norm_sampler(curr, stddev); // sample new candidate, place old value in params, @@ -82,8 +81,6 @@ inline void mh_posterior(ModelType& model, auto cand = norm_sampler(gen); params_it->next = curr; var.set_value(cand); - log_alpha -= dist::NormalBase::log_pdf(cand, curr, stddev); - log_alpha += dist::NormalBase::log_pdf(curr, cand, stddev); ++params_it; } @@ -93,45 +90,29 @@ inline void mh_posterior(ModelType& model, // compute next candidate log pdf and update log_alpha double cand_log_pdf = model.log_pdf(); log_alpha += cand_log_pdf; + bool accept = (std::log(unif_sampler(gen)) <= log_alpha); - // accept - if (std::log(unif_sampler(gen)) <= log_alpha) { - // "current" sample for next iteration is already in the variables - // simply append to storage - auto add_to_storage = [iter](auto& eq_node) { - auto& var = eq_node.get_variable(); - using var_t = std::decay_t; - using state_t = typename var_traits::state_t; - if (var.get_state() == state_t::parameter) { - auto storage = var.get_storage(); - storage[iter] = var.get_value(); - } - }; - model.traverse(add_to_storage); - - // update current log pdf for next iteration - curr_log_pdf = cand_log_pdf; - } - - // reject - else { - // "current" sample for next iteration must be moved back from - // params vector into variables. - // Then, simply append to storage - auto add_to_storage = [params_it, iter](auto& eq_node) mutable { - auto& var = eq_node.get_variable(); - using var_t = std::decay_t; - using state_t = typename var_traits::state_t; - if (var.get_state() == state_t::parameter) { + // If accept, "current" sample for next iteration is already in the variables + // so simply append to storage. + // Otherwise, "current" sample for next iteration must be moved back from + // params vector into variables. + auto add_to_storage = [params_it, iter, accept](auto& eq_node) mutable { + auto& var = eq_node.get_variable(); + using var_t = std::decay_t; + using state_t = typename util::var_traits::state_t; + if (var.get_state() == state_t::parameter) { + if (!accept) { var.set_value(params_it->next); - auto storage = var.get_storage(); - storage[iter] = var.get_value(); ++params_it; - } - }; - model.traverse(add_to_storage); - } + } + auto storage = var.get_storage(); + storage[iter] = var.get_value(); + } + }; + model.traverse(add_to_storage); + // update current log pdf for next iteration + if (accept) curr_log_pdf = cand_log_pdf; } } diff --git a/include/autoppl/distribution/dist_expr.hpp b/include/autoppl/distribution/dist_expr.hpp deleted file mode 100644 index 1eb23168..00000000 --- a/include/autoppl/distribution/dist_expr.hpp +++ /dev/null @@ -1,29 +0,0 @@ -#pragma once -#include - -namespace ppl { -namespace dist { - -template -struct DistExpr -{ - Derived& self() - { return static_cast(*this); } - - const Derived& self() const - { return static_cast(*this); } -}; - -template -inline constexpr bool is_dist_expr_v = - std::is_convertible_v>; - -#ifdef AUTOPPL_USE_CONCEPTS -// TODO: definition should be extended with a stronger -// restriction on T with interface checking. -template -concept dist_expressable = is_dist_expr_v; -#endif - -} // namespace dist -} // namespace ppl diff --git a/include/autoppl/expr_builder.hpp b/include/autoppl/expr_builder.hpp index 8fd0964f..53f50449 100644 --- a/include/autoppl/expr_builder.hpp +++ b/include/autoppl/expr_builder.hpp @@ -1,12 +1,12 @@ #pragma once -#include -#include -#include -#include -#include -#include -#include -#include +#include +#include +#include +#include +#include +#include +#include +#include namespace ppl { @@ -21,11 +21,51 @@ namespace ppl { namespace details { -using cont_raw_param_t = double; +/* + * Converter from arbitrary (decayed) type to valid continuous parameter type + * by the following mapping: + * - is_var_v true => VariableViewer + * - T is same as cont_raw_param_t => Constant + * - is_var_expr_v true => T + * Assumes each condition is non-overlapping. + */ +template +struct convert_to_cont_dist_param +{}; + +template +struct convert_to_cont_dist_param> && + !std::is_same_v, util::cont_raw_param_t> && + !util::is_var_expr_v> + >> +{ + using type = expr::VariableViewer>; +}; + +template +struct convert_to_cont_dist_param> && + std::is_same_v, util::cont_raw_param_t> && + !util::is_var_expr_v> + >> +{ + using type = expr::Constant>; +}; + +template +struct convert_to_cont_dist_param> && + !std::is_same_v, util::cont_raw_param_t> && + util::is_var_expr_v> + >> +{ + using type = T; +}; + template -inline constexpr bool is_cont_param_valid = - expr::is_var_expr_v || - std::is_convertible_v; +using convert_to_cont_dist_param_t = + typename convert_to_cont_dist_param::type; } // namespace details @@ -36,12 +76,16 @@ inline constexpr bool is_cont_param_valid = * See var_expr.hpp for more information. */ template -inline constexpr auto uniform(const MinType& min_expr, - const MaxType& max_expr) +inline constexpr auto uniform(MinType&& min_expr, + MaxType&& max_expr) { - static_assert(details::is_cont_param_valid); - static_assert(details::is_cont_param_valid); - return dist::Uniform(min_expr, max_expr); + using min_t = details::convert_to_cont_dist_param_t; + using max_t = details::convert_to_cont_dist_param_t; + + min_t wrap_min_expr = std::forward(min_expr); + max_t wrap_max_expr = std::forward(max_expr); + + return expr::Uniform(wrap_min_expr, wrap_max_expr); } #else #endif @@ -53,13 +97,18 @@ inline constexpr auto uniform(const MinType& min_expr, * See var_expr.hpp for more information. */ template -inline constexpr auto normal(const MeanType& mean_expr, - const StddevType& stddev_expr) +inline constexpr auto normal(MeanType&& mean_expr, + StddevType&& stddev_expr) { - static_assert(details::is_cont_param_valid); - static_assert(details::is_cont_param_valid); - return dist::Normal(mean_expr, stddev_expr); + using mean_t = details::convert_to_cont_dist_param_t; + using stddev_t = details::convert_to_cont_dist_param_t; + + mean_t wrap_mean_expr = std::forward(mean_expr); + stddev_t wrap_stddev_expr = std::forward(stddev_expr); + + return expr::Normal(wrap_mean_expr, wrap_stddev_expr); } + #else #endif @@ -72,8 +121,7 @@ inline constexpr auto normal(const MeanType& mean_expr, template inline constexpr auto bernoulli(const ProbType& p_expr) { - static_assert(details::is_cont_param_valid); - return dist::Bernoulli(p_expr); + return expr::Bernoulli(p_expr); } #else #endif @@ -90,11 +138,8 @@ inline constexpr auto bernoulli(const ProbType& p_expr) */ template inline constexpr auto operator|=(Variable& var, - const DistType& dist) -{ - static_assert(dist::is_dist_expr_v); - return expr::EqNode(var, dist); -} + DistType&& dist) +{ return expr::EqNode(var, std::forward(dist)); } #else #endif @@ -105,12 +150,11 @@ inline constexpr auto operator|=(Variable& var, * Ex. (x |= uniform(0,1), y |= uniform(0, 2)) */ template -inline constexpr auto operator,(const LHSNodeType& lhs, - const RHSNodeType& rhs) -{ - static_assert(expr::is_model_expr_v); - static_assert(expr::is_model_expr_v); - return expr::GlueNode(lhs, rhs); +inline constexpr auto operator,(LHSNodeType&& lhs, + RHSNodeType&& rhs) +{ + return expr::GlueNode(std::forward(lhs), + std::forward(rhs)); } #else #endif diff --git a/include/autoppl/distribution/bernoulli.hpp b/include/autoppl/expression/distribution/bernoulli.hpp similarity index 66% rename from include/autoppl/distribution/bernoulli.hpp rename to include/autoppl/expression/distribution/bernoulli.hpp index 18424b79..a6870f88 100644 --- a/include/autoppl/distribution/bernoulli.hpp +++ b/include/autoppl/expression/distribution/bernoulli.hpp @@ -1,20 +1,20 @@ #pragma once #include #include -#include -#include -#include -#include -#include +#include +#include +#include namespace ppl { -namespace dist { +namespace expr { template -struct Bernoulli : public DistExpr> +struct Bernoulli { - using value_t = uint8_t; - using param_value_t = typename var_traits::value_t; + static_assert(util::is_var_expr_v); + + using value_t = util::disc_raw_param_t; + using param_value_t = typename util::var_expr_traits::value_t; using dist_value_t = typename BernoulliBase::dist_value_t; Bernoulli(p_type p) @@ -39,5 +39,5 @@ struct Bernoulli : public DistExpr> p_type p_; }; -} // namespace dist +} // namespace expr } // namespace ppl diff --git a/include/autoppl/distribution/density.hpp b/include/autoppl/expression/distribution/density.hpp similarity index 97% rename from include/autoppl/distribution/density.hpp rename to include/autoppl/expression/distribution/density.hpp index 492d0827..ea267e1e 100644 --- a/include/autoppl/distribution/density.hpp +++ b/include/autoppl/expression/distribution/density.hpp @@ -1,8 +1,9 @@ #pragma once #include +#include namespace ppl { -namespace dist { +namespace expr { /* * The Base objects contain static member functions @@ -86,5 +87,5 @@ struct BernoulliBase } }; -} // namespace dist +} // namespace expr } // namespace ppl diff --git a/include/autoppl/distribution/normal.hpp b/include/autoppl/expression/distribution/normal.hpp similarity index 67% rename from include/autoppl/distribution/normal.hpp rename to include/autoppl/expression/distribution/normal.hpp index f069e693..b9970126 100644 --- a/include/autoppl/distribution/normal.hpp +++ b/include/autoppl/expression/distribution/normal.hpp @@ -1,22 +1,23 @@ #pragma once #include #include -#include -#include -#include -#include -#include +#include +#include +#include namespace ppl { -namespace dist { +namespace expr { template -struct Normal : public DistExpr> +struct Normal { - using value_t = double; + static_assert(util::is_var_expr_v); + static_assert(util::is_var_expr_v); + + using value_t = util::cont_raw_param_t; using param_value_t = std::common_type_t< - typename var_traits::value_t, - typename var_traits::value_t + typename util::var_expr_traits::value_t, + typename util::var_expr_traits::value_t >; using dist_value_t = typename NormalBase::dist_value_t; @@ -45,5 +46,5 @@ struct Normal : public DistExpr> stddev_type stddev_; }; -} // namespace dist +} // namespace expr } // namespace ppl diff --git a/include/autoppl/distribution/uniform.hpp b/include/autoppl/expression/distribution/uniform.hpp similarity index 68% rename from include/autoppl/distribution/uniform.hpp rename to include/autoppl/expression/distribution/uniform.hpp index 41c830c3..82b27bf2 100644 --- a/include/autoppl/distribution/uniform.hpp +++ b/include/autoppl/expression/distribution/uniform.hpp @@ -1,22 +1,23 @@ #pragma once #include #include -#include -#include -#include -#include -#include +#include +#include +#include namespace ppl { -namespace dist { +namespace expr { template -struct Uniform : public DistExpr> +struct Uniform { - using value_t = double; + static_assert(util::is_var_expr_v); + static_assert(util::is_var_expr_v); + + using value_t = util::cont_raw_param_t; using param_value_t = std::common_type_t< - typename var_traits::value_t, - typename var_traits::value_t + typename util::var_expr_traits::value_t, + typename util::var_expr_traits::value_t >; using dist_value_t = typename UniformBase::dist_value_t; @@ -45,5 +46,5 @@ struct Uniform : public DistExpr> max_type max_; }; -} // namespace dist +} // namespace expr } // namespace ppl diff --git a/include/autoppl/expression/model.hpp b/include/autoppl/expression/model/model.hpp similarity index 82% rename from include/autoppl/expression/model.hpp rename to include/autoppl/expression/model/model.hpp index 5e06c6ad..78adb7c0 100644 --- a/include/autoppl/expression/model.hpp +++ b/include/autoppl/expression/model/model.hpp @@ -2,33 +2,26 @@ #include #include #include -#include -#include +#include +#include +#include namespace ppl { namespace expr { -namespace details { - -template -struct IdentityVarFunctor -{ - using value_t = typename std::iterator_traits::value_type; - value_t& operator()(value_t& var) - { return var; } -}; - -} // namespace details /* * This class represents a "node" in the model expression * that relates a var with a distribution. */ template -struct EqNode : public ModelExpr> +struct EqNode { + static_assert(util::is_var_v); + static_assert(util::is_dist_expr_v); + using var_t = VarType; using dist_t = DistType; - using dist_value_t = typename dist_traits::dist_value_t; + using dist_value_t = typename util::dist_expr_traits::dist_value_t; EqNode(var_t& var, const dist_t& dist) noexcept @@ -77,13 +70,16 @@ struct EqNode : public ModelExpr> * "glues" two sub-model expressions. */ template -struct GlueNode : public ModelExpr> +struct GlueNode { + static_assert(util::is_model_expr_v); + static_assert(util::is_model_expr_v); + using left_node_t = LHSNodeType; using right_node_t = RHSNodeType; using dist_value_t = std::common_type_t< - typename node_traits::dist_value_t, - typename node_traits::dist_value_t + typename util::model_expr_traits::dist_value_t, + typename util::model_expr_traits::dist_value_t >; GlueNode(const left_node_t& lhs, diff --git a/include/autoppl/expression/model_expr.hpp b/include/autoppl/expression/model_expr.hpp deleted file mode 100644 index 05a26304..00000000 --- a/include/autoppl/expression/model_expr.hpp +++ /dev/null @@ -1,29 +0,0 @@ -#pragma once -#include - -namespace ppl { -namespace expr { - -template -struct ModelExpr -{ - Derived& self() - { return static_cast(*this); } - - const Derived& self() const - { return static_cast(*this); } -}; - -template -inline constexpr bool is_model_expr_v = - std::is_convertible_v>; - -#ifdef AUTOPPL_USE_CONCEPTS -// TODO: definition should be extended with a stronger -// restriction on T with interface checking. -template -concept model_expressable = is_model_expr_v; -#endif - -} // namespace expr -} // namespace ppl diff --git a/include/autoppl/expression/var_expr.hpp b/include/autoppl/expression/var_expr.hpp deleted file mode 100644 index c891e5b9..00000000 --- a/include/autoppl/expression/var_expr.hpp +++ /dev/null @@ -1,29 +0,0 @@ -#pragma once -#include - -namespace ppl { -namespace expr { - -template -struct VarExpr -{ - Derived& self() - { return static_cast(*this); } - - const Derived& self() const - { return static_cast(*this); } -}; - -template -inline constexpr bool is_var_expr_v = - std::is_convertible_v>; - -#ifdef AUTOPPL_USE_CONCEPTS -// TODO: definition should be extended with a stronger -// restriction on T with interface checking. -template -concept var_expressable = is_var_expr_v; -#endif - -} // namespace expr -} // namespace ppl diff --git a/include/autoppl/expression/variable/constant.hpp b/include/autoppl/expression/variable/constant.hpp new file mode 100644 index 00000000..28c95856 --- /dev/null +++ b/include/autoppl/expression/variable/constant.hpp @@ -0,0 +1,20 @@ +#pragma once + +namespace ppl { +namespace expr { + +template +struct Constant +{ + using value_t = ValueType; + Constant(value_t c) + : c_{c} + {} + operator value_t() const { return c_; } + +private: + value_t c_; +}; + +} // namespace expr +} // namespace ppl diff --git a/include/autoppl/expression/variable/variable_viewer.hpp b/include/autoppl/expression/variable/variable_viewer.hpp new file mode 100644 index 00000000..9e137913 --- /dev/null +++ b/include/autoppl/expression/variable/variable_viewer.hpp @@ -0,0 +1,33 @@ +#pragma once +#include +#include + +namespace ppl { +namespace expr { + +/* + * VariableViewer is a viewer of some variable type. + * It will mainly be used to view Variable class defined in autoppl/variable.hpp. + */ +template +struct VariableViewer +{ + static_assert(util::is_var_v); + + using var_t = VariableType; + using value_t = typename util::var_traits::value_t; + + VariableViewer(var_t& var) + : var_ref_{var} + {} + + operator value_t() const + { return static_cast(var_ref_.get()); } + +private: + using var_ref_t = std::reference_wrapper; + var_ref_t var_ref_; +}; + +} // namespace expr +} // namespace ppl diff --git a/include/autoppl/util/concept.hpp b/include/autoppl/util/concept.hpp new file mode 100644 index 00000000..8ccaafe7 --- /dev/null +++ b/include/autoppl/util/concept.hpp @@ -0,0 +1,127 @@ +#pragma once +#include + +/* + * Metaprogramming tool to check if name is a (public) + * member alias of a given type T. + * All versions must be placed in this file for ease of maintenance. + * Macro definition is undefined at the end of the file. + * + * Ex. with "name" as "value_t" + * + namespace details { + template + struct has_type_value_t + { + private: + template static void impl(decltype(typename V::value_t(), int())); + template static bool impl(char); + public: + static constexpr bool value = std::is_same(0))>::value; + }; + } + template + inline constexpr bool has_type_value_t_v = + details::has_type_value_t::value; + */ + +#define DEFINE_HAS_TYPE(name) \ + namespace details { \ + template \ + struct has_type_##name \ + { \ + private: \ + template static void impl(typename V::name*); \ + template static bool impl(...); \ + public: \ + static constexpr bool value = std::is_same(0))>::value; \ + }; \ + \ + template \ + struct get_type_##name \ + { \ + using type = invalid_tag; \ + }; \ + template \ + struct get_type_##name \ + { \ + using type = typename T::name; \ + }; \ + } \ + template \ + inline constexpr bool has_type_##name##_v = \ + details::has_type_##name::value; \ + template \ + using get_type_##name##_t = \ + typename details::get_type_##name>::type; + +/* + * Metaprogramming tool to check if name is a (public) + * member function of a given type T. + * All versions must be placed in this file for ease of maintenance. + * Macro definition is undefined at the end of the file. + * + * Ex. with "name" as "pdf" + * + namespace details { + template + struct has_func_pdf + { + private: + template static void impl(decltype(&V::pdf)); + template static bool impl(...); + public: + static constexpr bool value = std::is_same(0))>::value; + }; + } + template + inline constexpr bool has_func_pdf_v = + details::has_func_pdf::value; + */ + +#define DEFINE_HAS_FUNC(name) \ + namespace details { \ + template \ + struct has_func_##name \ + { \ + private: \ + template static void impl(decltype(&V::name)); \ + template static bool impl(...); \ + public: \ + static constexpr bool value = std::is_same(0))>::value; \ + }; \ + } \ + template \ + inline constexpr bool has_func_##name##_v = \ + details::has_func_##name::value; + +namespace ppl { +namespace util { + +struct invalid_tag {}; + +DEFINE_HAS_TYPE(value_t); +DEFINE_HAS_TYPE(pointer_t); +DEFINE_HAS_TYPE(const_pointer_t); +DEFINE_HAS_TYPE(state_t); + +DEFINE_HAS_TYPE(dist_value_t); + +DEFINE_HAS_FUNC(set_value); +DEFINE_HAS_FUNC(get_value); +DEFINE_HAS_FUNC(set_storage); +DEFINE_HAS_FUNC(get_storage); +DEFINE_HAS_FUNC(set_state); +DEFINE_HAS_FUNC(get_state); + +DEFINE_HAS_FUNC(pdf); +DEFINE_HAS_FUNC(log_pdf); + +DEFINE_HAS_FUNC(get_variable); +DEFINE_HAS_FUNC(get_distribution); + +} // namespace util +} // namespace ppl + +#undef DEFINE_HAS_FUNC +#undef DEFINE_HAS_TYPE diff --git a/include/autoppl/util/dist_expr_traits.hpp b/include/autoppl/util/dist_expr_traits.hpp new file mode 100644 index 00000000..e674e848 --- /dev/null +++ b/include/autoppl/util/dist_expr_traits.hpp @@ -0,0 +1,60 @@ +#pragma once +#include +#include + +namespace ppl { +namespace util { + +/* + * TODO: Samplable distribution expression concept? + */ + +/* + * TODO: continuous/discrete distribution expression concept? + */ + +/* + * Continuous distribution expressions can be constructed with this type. + */ +using cont_raw_param_t = double; + +/* + * Discrete distribution expressions can be constructed with this type. + */ +using disc_raw_param_t = int64_t; + +/* + * Traits for Distribution Expression classes. + * value_t type of value Variable represents during computation + * dist_value_t type of pdf/log_pdf value + */ +template +struct dist_expr_traits +{ + using value_t = typename DistExprType::value_t; + using dist_value_t = typename DistExprType::dist_value_t; +}; + +/* + * A distribution expression is any class that the following: + * - dist_expr_traits must be well-defined for T + * - T must have member function pdf + * - T must have member function log_pdf + */ +template +inline constexpr bool is_dist_expr_v = + has_type_value_t_v && + has_type_dist_value_t_v && + has_func_pdf_v && + has_func_log_pdf_v + ; + +#ifdef AUTOPPL_USE_CONCEPTS +// TODO: definition should be extended with a stronger +// restriction on T with interface checking. +template +concept dist_expressable = is_dist_expr_v; +#endif + +} // namespace util +} // namespace ppl diff --git a/include/autoppl/util/model_expr_traits.hpp b/include/autoppl/util/model_expr_traits.hpp new file mode 100644 index 00000000..f631b4cf --- /dev/null +++ b/include/autoppl/util/model_expr_traits.hpp @@ -0,0 +1,43 @@ +#pragma once +#include + +namespace ppl { +namespace util { + +/* + * Traits for Model Expression classes. + * dist_value_t type of value Variable represents during computation + */ +template +struct model_expr_traits +{ + using dist_value_t = typename NodeType::dist_value_t; +}; + +// TODO: +// - pdf and log_pdf remove from interface? +// - how to check if template member function exists (for traverse)? +template +inline constexpr bool is_model_expr_v = + has_type_dist_value_t_v && + has_func_pdf_v && + has_func_log_pdf_v + ; + +// TODO: not used currently +template +inline constexpr bool is_eq_node_expr_v = + is_model_expr_v && + has_func_get_variable_v && + has_func_get_distribution_v + ; + +#ifdef AUTOPPL_USE_CONCEPTS +// TODO: definition should be extended with a stronger +// restriction on T with interface checking. +template +concept model_expressable = is_model_expr_v; +#endif + +} // namespace util +} // namespace ppl diff --git a/include/autoppl/util/traits.hpp b/include/autoppl/util/traits.hpp index 5e812978..9c137cf2 100644 --- a/include/autoppl/util/traits.hpp +++ b/include/autoppl/util/traits.hpp @@ -1,7 +1,4 @@ #pragma once -#include - -namespace ppl { /* * The following classes list member aliases that @@ -9,36 +6,7 @@ namespace ppl { * Users should rely on these classes to grab member aliases. */ -template -struct var_traits -{ - using value_t = typename VarType::value_t; - using pointer_t = typename VarType::pointer_t; - using state_t = typename VarType::state_t; - - // TODO may have to move this to a different class for compile-time checking - static_assert(std::is_convertible_v); -}; - -// Specialization: when double or int, considered "trivial" variable. -// TODO: this was a quick fix for generalizing distribution value_t. -template <> -struct var_traits -{ - using value_t = double; -}; - -template -struct dist_traits -{ - using value_t = typename DistType::value_t; // generated value type - using dist_value_t = typename DistType::dist_value_t; // pdf value type -}; - -template -struct node_traits -{ - using dist_value_t = typename NodeType::dist_value_t; -}; - -} // namespace ppl +#include +#include +#include +#include diff --git a/include/autoppl/util/var_expr_traits.hpp b/include/autoppl/util/var_expr_traits.hpp new file mode 100644 index 00000000..8513448e --- /dev/null +++ b/include/autoppl/util/var_expr_traits.hpp @@ -0,0 +1,47 @@ +#pragma once +#include +#include + +namespace ppl { +namespace util { + +/* + * Traits for Variable Expression classes. + * value_t type of value Variable represents during computation + */ +template +struct var_expr_traits +{ + using value_t = typename VarExprType::value_t; +}; + +// Specialization: when double or int, considered "trivial" variable. +// TODO: this was a quick fix for generalizing distribution value_t. +template <> +struct var_expr_traits +{ + using value_t = double; +}; + +/* + * A variable expression is any class that the following: + * - is_var_v must be false + * - var_expr_traits must be well-defined for T + * - T must be convertible to its value_t + */ +template +inline constexpr bool is_var_expr_v = + !is_var_v && + has_type_value_t_v && + std::is_convertible_v> + ; + +#ifdef AUTOPPL_USE_CONCEPTS +// TODO: definition should be extended with a stronger +// restriction on T with interface checking. +template +concept var_expressable = is_var_expr_v; +#endif + +} // namespace util +} // namespace ppl diff --git a/include/autoppl/util/var_traits.hpp b/include/autoppl/util/var_traits.hpp new file mode 100644 index 00000000..173d0818 --- /dev/null +++ b/include/autoppl/util/var_traits.hpp @@ -0,0 +1,43 @@ +#pragma once +#include +#include + +namespace ppl { +namespace util { + +/* + * Traits for Variable-like classes. + * value_t type of value Variable represents during computation + * pointer_t storage pointer type + * state_t type of enum class state; must have "data" and "parameter" + */ +template +struct var_traits +{ + using value_t = typename VarType::value_t; + using pointer_t = typename VarType::pointer_t; + using state_t = typename VarType::state_t; +}; + +/* + * C++17 version of concepts to check var properties. + * - var_traits must be well-defined under type T + * - T must be convertible to its value_t + * - not possible to get overloads + */ +template +inline constexpr bool is_var_v = + has_type_value_t_v && + has_type_pointer_t_v && + has_type_const_pointer_t_v && + has_type_state_t_v && + has_func_set_value_v && + has_func_get_value_v && + has_func_set_storage_v && + has_func_set_state_v && + has_func_get_state_v && + std::is_convertible_v> + ; + +} // namespace util +} // namespace ppl diff --git a/include/autoppl/expression/variable.hpp b/include/autoppl/variable.hpp similarity index 90% rename from include/autoppl/expression/variable.hpp rename to include/autoppl/variable.hpp index eeca0092..1217124c 100644 --- a/include/autoppl/expression/variable.hpp +++ b/include/autoppl/variable.hpp @@ -1,5 +1,6 @@ #pragma once -#include +#include +#include namespace ppl { @@ -19,7 +20,7 @@ enum class var_state : bool { * a model expression and the users, who must supply storage of values associated with this var. */ template -struct Variable : public expr::VarExpr> +struct Variable { using value_t = ValueType; using pointer_t = value_t*; @@ -77,7 +78,7 @@ struct Variable : public expr::VarExpr> }; // Useful aliases -using cont_var = Variable; // continuous RV var -using disc_var = Variable; // discrete RV var +using cont_var = Variable; // continuous RV var +using disc_var = Variable; // discrete RV var } // namespace ppl diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 0383af60..a898ecdf 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -6,36 +6,83 @@ if (AUTOPPL_ENABLE_TEST_COVERAGE) endif() ###################################################### -# Expression Test +# Util Test ###################################################### -add_executable(expression_unittest - ${CMAKE_CURRENT_SOURCE_DIR}/expression/model_unittest.cpp +add_executable(util_unittest + ${CMAKE_CURRENT_SOURCE_DIR}/util/concept_unittest.cpp + ${CMAKE_CURRENT_SOURCE_DIR}/util/dist_expr_traits_unittest.cpp + ${CMAKE_CURRENT_SOURCE_DIR}/util/var_expr_traits_unittest.cpp + ${CMAKE_CURRENT_SOURCE_DIR}/util/var_traits_unittest.cpp + ) +target_compile_options(util_unittest PRIVATE -g -Wall -Werror -Wextra) +target_include_directories(util_unittest PRIVATE + ${GTEST_DIR}/include + ${CMAKE_CURRENT_SOURCE_DIR} + ) +if (AUTOPPL_ENABLE_TEST_COVERAGE) + target_link_libraries(util_unittest gcov) +endif() +target_link_libraries(util_unittest gtest_main pthread ${PROJECT_NAME}) +add_test(util_unittest util_unittest) + +###################################################### +# Variable Expression Test +###################################################### + +add_executable(var_expr_unittest + ${CMAKE_CURRENT_SOURCE_DIR}/expression/variable/variable_viewer_unittest.cpp + ${CMAKE_CURRENT_SOURCE_DIR}/expression/variable/constant_unittest.cpp + ) +target_compile_options(var_expr_unittest PRIVATE -g -Wall -Werror -Wextra) +target_include_directories(var_expr_unittest PRIVATE + ${GTEST_DIR}/include + ${CMAKE_CURRENT_SOURCE_DIR} ) -target_compile_options(expression_unittest PRIVATE -g -Wall -Werror -Wextra) -target_include_directories(expression_unittest PRIVATE ${GTEST_DIR}/include) if (AUTOPPL_ENABLE_TEST_COVERAGE) - target_link_libraries(expression_unittest gcov) + target_link_libraries(var_expr_unittest gcov) endif() -target_link_libraries(expression_unittest gtest_main pthread ${PROJECT_NAME}) -add_test(expression_unittest expression_unittest) +target_link_libraries(var_expr_unittest gtest_main pthread ${PROJECT_NAME}) +add_test(var_expr_unittest var_expr_unittest) ###################################################### -# Distribution Test +# Distribution Expression Test ###################################################### -add_executable(distribution_unittest - ${CMAKE_CURRENT_SOURCE_DIR}/distribution/uniform_unittest.cpp - ${CMAKE_CURRENT_SOURCE_DIR}/distribution/normal_unittest.cpp - ${CMAKE_CURRENT_SOURCE_DIR}/distribution/bernoulli_unittest.cpp +add_executable(dist_expr_unittest + ${CMAKE_CURRENT_SOURCE_DIR}/expression/distribution/bernoulli_unittest.cpp + ${CMAKE_CURRENT_SOURCE_DIR}/expression/distribution/density_unittest.cpp + ${CMAKE_CURRENT_SOURCE_DIR}/expression/distribution/normal_unittest.cpp + ${CMAKE_CURRENT_SOURCE_DIR}/expression/distribution/uniform_unittest.cpp + ) +target_compile_options(dist_expr_unittest PRIVATE -g -Wall -Werror -Wextra) +target_include_directories(dist_expr_unittest PRIVATE + ${GTEST_DIR}/include + ${CMAKE_CURRENT_SOURCE_DIR} ) -target_compile_options(distribution_unittest PRIVATE -g -Wall -Werror -Wextra) -target_include_directories(distribution_unittest PRIVATE ${GTEST_DIR}/include) if (AUTOPPL_ENABLE_TEST_COVERAGE) - target_link_libraries(distribution_unittest gcov) + target_link_libraries(dist_expr_unittest gcov) endif() -target_link_libraries(distribution_unittest gtest_main pthread ${PROJECT_NAME}) -add_test(distribution_unittest distribution_unittest) +target_link_libraries(dist_expr_unittest gtest_main pthread ${PROJECT_NAME}) +add_test(dist_expr_unittest dist_expr_unittest) + +###################################################### +# Model Expression Test +###################################################### + +add_executable(model_expr_unittest + ${CMAKE_CURRENT_SOURCE_DIR}/expression/model/model_unittest.cpp + ) +target_compile_options(model_expr_unittest PRIVATE -g -Wall -Werror -Wextra) +target_include_directories(model_expr_unittest PRIVATE + ${GTEST_DIR}/include + ${CMAKE_CURRENT_SOURCE_DIR} + ) +if (AUTOPPL_ENABLE_TEST_COVERAGE) + target_link_libraries(model_expr_unittest gcov) +endif() +target_link_libraries(model_expr_unittest gtest_main pthread ${PROJECT_NAME}) +add_test(model_expr_unittest model_expr_unittest) ###################################################### # Algorithm Test @@ -45,9 +92,30 @@ add_executable(algorithm_unittest ${CMAKE_CURRENT_SOURCE_DIR}/algorithm/mh_unittest.cpp ) target_compile_options(algorithm_unittest PRIVATE -g -Wall -Werror -Wextra) -target_include_directories(algorithm_unittest PRIVATE ${GTEST_DIR}/include) +target_include_directories(algorithm_unittest PRIVATE + ${GTEST_DIR}/include + ${CMAKE_CURRENT_SOURCE_DIR} + ) if (AUTOPPL_ENABLE_TEST_COVERAGE) target_link_libraries(algorithm_unittest gcov) endif() target_link_libraries(algorithm_unittest gtest_main pthread ${PROJECT_NAME}) add_test(algorithm_unittest algorithm_unittest) + +###################################################### +# Expression Builder Test +###################################################### + +add_executable(expr_builder_unittest + ${CMAKE_CURRENT_SOURCE_DIR}/expr_builder_unittest.cpp + ) +target_compile_options(expr_builder_unittest PRIVATE -g -Wall -Werror -Wextra) +target_include_directories(expr_builder_unittest PRIVATE + ${GTEST_DIR}/include + ${CMAKE_CURRENT_SOURCE_DIR} + ) +if (AUTOPPL_ENABLE_TEST_COVERAGE) + target_link_libraries(expr_builder_unittest gcov) +endif() +target_link_libraries(expr_builder_unittest gtest_main pthread ${PROJECT_NAME}) +add_test(expr_builder_unittest expr_builder_unittest) diff --git a/test/algorithm/mh_unittest.cpp b/test/algorithm/mh_unittest.cpp index c1299d10..b88400ac 100644 --- a/test/algorithm/mh_unittest.cpp +++ b/test/algorithm/mh_unittest.cpp @@ -1,41 +1,15 @@ #include "gtest/gtest.h" -#include #include +#include #include #include +#include namespace ppl { -template -inline void plot_hist(const ArrayType& arr, - double step_size = .5) -{ - constexpr size_t nstars = 100; // maximum number of stars to distribute - - const int64_t min = std::floor(*std::min_element(arr.begin(), arr.end())); - const int64_t max = std::ceil(*std::max_element(arr.begin(), arr.end())); - const uint64_t range = max - min; - const uint64_t n_hist = std::floor(range/step_size); - - // keeps count for each histogram bar - std::vector counter(n_hist, 0); - - for (auto x : arr) { - ++counter[std::floor((x - min) / step_size)]; - } - - EXPECT_EQ(std::accumulate(counter.begin(), counter.end(), 0), (int) arr.size()); - - for (size_t i = 0; i < n_hist; ++i) { - std::cout << i << "-" << (i+1) << ": " << '\t'; - std::cout << std::string(counter[i] * nstars/arr.size(), '*') << std::endl; - } - -} - TEST(mh_unittest, plot_hist_sanity) { - static constexpr size_t sample_size = 5000; + static constexpr size_t sample_size = 20000; std::array storage = {0.}; std::normal_distribution normal_sampler(0., 1.); std::mt19937 gen; @@ -45,16 +19,28 @@ TEST(mh_unittest, plot_hist_sanity) plot_hist(storage); } +/* + * Fixture for Metropolis-Hastings + */ struct mh_fixture : ::testing::Test { protected: - static constexpr size_t sample_size = 5000; + static constexpr size_t sample_size = 20000; std::array storage = {0.}; - Variable theta; + Variable theta, x; mh_fixture() : theta{storage.data()} {} + + double sample_average(size_t burn) + { + double sum = std::accumulate( + std::next(storage.begin(), burn), + storage.end(), + 0.); + return sum / (storage.size() - burn); + } }; TEST_F(mh_fixture, sample_std_normal) @@ -62,13 +48,27 @@ TEST_F(mh_fixture, sample_std_normal) auto model = (theta |= normal(0., 1.)); mh_posterior(model, sample_size); plot_hist(storage); + EXPECT_NEAR(sample_average(1000), 0., 0.01); } TEST_F(mh_fixture, sample_uniform) { auto model = (theta |= uniform(0., 1.)); mh_posterior(model, sample_size); - plot_hist(storage, 0.05); + plot_hist(storage, 0.05, 0., 1.); + EXPECT_NEAR(sample_average(1000), 0.5, 0.01); +} + +TEST_F(mh_fixture, sample_unif_normal_posterior) +{ + x.observe(3.); + auto model = ( + theta |= uniform(-20., 20.), + x |= normal(theta, 1.) + ); + mh_posterior(model, sample_size); + plot_hist(storage); + EXPECT_NEAR(sample_average(1000), 3.0, 0.01); } } // namespace ppl diff --git a/test/distribution/bernoulli_unittest.cpp b/test/distribution/bernoulli_unittest.cpp deleted file mode 100644 index 267ef18e..00000000 --- a/test/distribution/bernoulli_unittest.cpp +++ /dev/null @@ -1,44 +0,0 @@ -#include -#include - -#include -#include - -#include "gtest/gtest.h" - -namespace ppl { -namespace dist { - -struct bernoulli_dist_fixture : ::testing::Test { -protected: - Variable x {0.6}; - - Bernoulli dist1 = Bernoulli(0.6); - Bernoulli > dist2 = Bernoulli(x); -}; - -TEST_F(bernoulli_dist_fixture, sanity_bernoulli_test) { - EXPECT_EQ(dist1.p(), 0.6); - EXPECT_EQ(dist2.p(), 0.6); -} - -TEST_F(bernoulli_dist_fixture, simple_bernoulli) { - EXPECT_DOUBLE_EQ(dist1.pdf(1), dist1.p()); - EXPECT_DOUBLE_EQ(dist1.pdf(1), 0.6); - EXPECT_DOUBLE_EQ(dist1.pdf(0), 1 - dist1.p()); - EXPECT_DOUBLE_EQ(dist1.pdf(0), 0.4); - EXPECT_DOUBLE_EQ(dist1.pdf(2), 0.0); -} - -TEST_F(bernoulli_dist_fixture, bernoulli_sampling) { - std::random_device rd{}; - std::mt19937 gen{rd()}; - - for (int i = 0; i < 100; i++) { - int sample = dist1.sample(gen); - EXPECT_TRUE(sample == 0 || sample == 1); - } -} - -} // namespace dist -} // namespace ppl diff --git a/test/distribution/normal_unittest.cpp b/test/distribution/normal_unittest.cpp deleted file mode 100644 index 3dea0eff..00000000 --- a/test/distribution/normal_unittest.cpp +++ /dev/null @@ -1,48 +0,0 @@ -#include -#include - -#include -#include - -#include "gtest/gtest.h" - -namespace ppl { -namespace dist { - -struct normal_dist_fixture : ::testing::Test { -protected: - Variable mu {0.}; - Variable sigma {1.}; - Normal dist1 = Normal(0., 1.); - Normal, Variable > dist2 = Normal(mu, sigma); -}; - -TEST_F(normal_dist_fixture, sanity_normal_test) { - EXPECT_EQ(dist1.mean(), 0.0); - EXPECT_EQ(dist1.stddev(), 1.0); - - EXPECT_EQ(dist2.mean(), 0.0); - EXPECT_EQ(dist2.stddev(), 1.0); -} - -TEST_F(normal_dist_fixture, simple_gaussian) { - EXPECT_DOUBLE_EQ(dist1.pdf(0.0), 0.3989422804014327); - EXPECT_DOUBLE_EQ(dist1.pdf(-0.5), 0.3520653267642995); - EXPECT_DOUBLE_EQ(dist1.pdf(4), 0.00013383022576488537); - - EXPECT_DOUBLE_EQ(dist1.log_pdf(0.0), std::log(dist1.pdf(0.0))); - EXPECT_DOUBLE_EQ(dist1.log_pdf(-0.5), std::log(dist1.pdf(-0.5))); - EXPECT_DOUBLE_EQ(dist1.log_pdf(4), std::log(dist1.pdf(4))); - - - EXPECT_DOUBLE_EQ(dist2.pdf(0.0), 0.3989422804014327); - EXPECT_DOUBLE_EQ(dist2.pdf(-0.5), 0.3520653267642995); - EXPECT_DOUBLE_EQ(dist2.pdf(4), 0.00013383022576488537); - - EXPECT_DOUBLE_EQ(dist2.log_pdf(0.0), std::log(dist2.pdf(0.0))); - EXPECT_DOUBLE_EQ(dist2.log_pdf(-0.5), std::log(dist2.pdf(-0.5))); - EXPECT_DOUBLE_EQ(dist2.log_pdf(4), std::log(dist1.pdf(4))); -} - -} // namespace dist -} // namespace ppl diff --git a/test/distribution/uniform_unittest.cpp b/test/distribution/uniform_unittest.cpp deleted file mode 100644 index cf1d7b4b..00000000 --- a/test/distribution/uniform_unittest.cpp +++ /dev/null @@ -1,54 +0,0 @@ -#include -#include - -#include -#include - -#include "gtest/gtest.h" - -namespace ppl { -namespace dist { - -struct uniform_dist_fixture : ::testing::Test { -protected: - Variable x {0.5}; - Variable y {0.1}; - Uniform dist1 = Uniform(0., 1.); - Uniform > dist2 = Uniform(0., x); - Uniform, Variable > dist3 = Uniform(y, x); -}; - -TEST_F(uniform_dist_fixture, sanity_uniform_test) { - EXPECT_EQ(dist1.min(), 0.0); - EXPECT_EQ(dist1.max(), 1.0); - - EXPECT_EQ(dist2.min(), 0.0); - EXPECT_EQ(dist2.max(), 0.5); - - EXPECT_EQ(dist3.min(), 0.1); - EXPECT_EQ(dist3.max(), 0.5); -} - -TEST_F(uniform_dist_fixture, simple_uniform) { - EXPECT_DOUBLE_EQ(dist1.pdf(1.1), 0.0); - EXPECT_DOUBLE_EQ(dist1.pdf(1.0), 0.0); - - EXPECT_DOUBLE_EQ(dist2.pdf(0.25), 2.0); - EXPECT_DOUBLE_EQ(dist2.pdf(-0.1), 0.0); - - EXPECT_DOUBLE_EQ(dist3.pdf(0.25), 2.5); -} - -TEST_F(uniform_dist_fixture, uniform_sampling) { - std::random_device rd{}; - std::mt19937 gen{rd()}; - - for (int i = 0; i < 100; i++) { - double sample = dist1.sample(gen); - EXPECT_GT(sample, 0.0); - EXPECT_LT(sample, 1.0); - } -} - -} // namespace dist -} // namespace ppl diff --git a/test/expr_builder_unittest.cpp b/test/expr_builder_unittest.cpp new file mode 100644 index 00000000..5e868895 --- /dev/null +++ b/test/expr_builder_unittest.cpp @@ -0,0 +1,55 @@ +#include "gtest/gtest.h" +#include +#include + +namespace ppl { + +struct expr_builder_fixture : ::testing::Test +{ +protected: +}; + +TEST_F(expr_builder_fixture, convert_to_cont_dist_param_var) +{ + using namespace details; + static_assert(std::is_same_v>); + static_assert(util::is_var_v); + static_assert(!std::is_same_v); + static_assert(!util::is_var_expr_v); + static_assert(std::is_same_v< + convert_to_cont_dist_param_t, + expr::VariableViewer + >); +} + +TEST_F(expr_builder_fixture, convert_to_cont_dist_param_raw) +{ + using namespace details; + using data_t = util::cont_raw_param_t; + static_assert(std::is_same_v>); + static_assert(!util::is_var_v); + static_assert(std::is_same_v); + static_assert(!util::is_var_expr_v); + static_assert(std::is_same_v< + convert_to_cont_dist_param_t, + expr::Constant + >); +} + +TEST_F(expr_builder_fixture, convert_to_cont_dist_param_var_expr) +{ + using namespace details; + static_assert(!util::is_var_v); + static_assert(!std::is_same_v); + static_assert(util::is_var_expr_v); + static_assert(std::is_same_v< + convert_to_cont_dist_param_t, + MockVarExpr& + >); + static_assert(std::is_same_v< + convert_to_cont_dist_param_t, + MockVarExpr&& + >); +} + +} // namespace ppl diff --git a/test/expression/distribution/bernoulli_unittest.cpp b/test/expression/distribution/bernoulli_unittest.cpp new file mode 100644 index 00000000..19641b4c --- /dev/null +++ b/test/expression/distribution/bernoulli_unittest.cpp @@ -0,0 +1,65 @@ +#include "gtest/gtest.h" +#include +#include +#include +#include +#include + +namespace ppl { +namespace expr { + +struct bernoulli_fixture : ::testing::Test +{ +protected: + using value_t = typename MockVarExpr::value_t; + static constexpr size_t sample_size = 1000; + double p = 0.6; + MockVarExpr x{p}; + Bernoulli bern = {x}; + std::array sample = {0.}; +}; + +TEST_F(bernoulli_fixture, ctor) +{ + static_assert(util::is_dist_expr_v>); +} + +TEST_F(bernoulli_fixture, bernoulli_check_params) { + EXPECT_DOUBLE_EQ(bern.p(), static_cast(x)); +} + +TEST_F(bernoulli_fixture, bernoulli_pdf_delegate) { + EXPECT_DOUBLE_EQ(bern.pdf(-10), BernoulliBase::pdf(-10, p)); + EXPECT_DOUBLE_EQ(bern.pdf(-7), BernoulliBase::pdf(-7, p)); + EXPECT_DOUBLE_EQ(bern.pdf(-3), BernoulliBase::pdf(-3, p)); + EXPECT_DOUBLE_EQ(bern.pdf(0), BernoulliBase::pdf(0, p)); + EXPECT_DOUBLE_EQ(bern.pdf(1), BernoulliBase::pdf(1, p)); + EXPECT_DOUBLE_EQ(bern.pdf(3), BernoulliBase::pdf(3, p)); + EXPECT_DOUBLE_EQ(bern.pdf(6), BernoulliBase::pdf(6, p)); + EXPECT_DOUBLE_EQ(bern.pdf(16), BernoulliBase::pdf(16, p)); +} + +TEST_F(bernoulli_fixture, bernoulli_log_pdf_delegate) { + EXPECT_DOUBLE_EQ(bern.log_pdf(-10), BernoulliBase::log_pdf(-10, p)); + EXPECT_DOUBLE_EQ(bern.log_pdf(-7), BernoulliBase::log_pdf(-7, p)); + EXPECT_DOUBLE_EQ(bern.log_pdf(-3), BernoulliBase::log_pdf(-3, p)); + EXPECT_DOUBLE_EQ(bern.log_pdf(0), BernoulliBase::log_pdf(0, p)); + EXPECT_DOUBLE_EQ(bern.log_pdf(1), BernoulliBase::log_pdf(1, p)); + EXPECT_DOUBLE_EQ(bern.log_pdf(3), BernoulliBase::log_pdf(3, p)); + EXPECT_DOUBLE_EQ(bern.log_pdf(6), BernoulliBase::log_pdf(6, p)); + EXPECT_DOUBLE_EQ(bern.log_pdf(16), BernoulliBase::log_pdf(16, p)); +} + +TEST_F(bernoulli_fixture, bernoulli_sample) { + std::random_device rd{}; + std::mt19937 gen{rd()}; + + for (size_t i = 0; i < sample_size; i++) { + sample[i] = bern.sample(gen); + } + + plot_hist(sample); +} + +} // namespace expr +} // namespace ppl diff --git a/test/expression/distribution/density_unittest.cpp b/test/expression/distribution/density_unittest.cpp new file mode 100644 index 00000000..6af92576 --- /dev/null +++ b/test/expression/distribution/density_unittest.cpp @@ -0,0 +1,157 @@ +#include "gtest/gtest.h" +#include +#include + +namespace ppl { +namespace expr { + +struct density_fixture : ::testing::Test +{ +protected: + double x = 0.; + double min = -2.3; + double max = 2.7; + double mean = 0.3; + double stddev = 1.3; + double tol = 1e-15; + double p = 0.41; +}; + +/* + * Continuous distribution + */ + +TEST_F(density_fixture, uniform_pdf_in_range) +{ + EXPECT_DOUBLE_EQ(UniformBase::pdf(-2.2999999999, min, max), 0.2); + EXPECT_DOUBLE_EQ(UniformBase::pdf(-2., min, max), 0.2); + EXPECT_DOUBLE_EQ(UniformBase::pdf(-1.423, min, max), 0.2); + EXPECT_DOUBLE_EQ(UniformBase::pdf(0., min, max), 0.2); + EXPECT_DOUBLE_EQ(UniformBase::pdf(1.31, min, max), 0.2); + EXPECT_DOUBLE_EQ(UniformBase::pdf(2.41, min, max), 0.2); + EXPECT_DOUBLE_EQ(UniformBase::pdf(2.69999999999, min, max), 0.2); +} + +TEST_F(density_fixture, uniform_pdf_out_of_range) +{ + EXPECT_DOUBLE_EQ(UniformBase::pdf(-100, min, max), 0.); + EXPECT_DOUBLE_EQ(UniformBase::pdf(-3.41, min, max), 0.); + EXPECT_DOUBLE_EQ(UniformBase::pdf(-2.3, min, max), 0.); + EXPECT_DOUBLE_EQ(UniformBase::pdf(2.7, min, max), 0.); + EXPECT_DOUBLE_EQ(UniformBase::pdf(3.5, min, max), 0.); + EXPECT_DOUBLE_EQ(UniformBase::pdf(3214, min, max), 0.); +} + +TEST_F(density_fixture, uniform_log_pdf_in_range) +{ + EXPECT_DOUBLE_EQ(UniformBase::log_pdf(-2.2999999999, min, max), std::log(0.2)); + EXPECT_DOUBLE_EQ(UniformBase::log_pdf(-2., min, max), std::log(0.2)); + EXPECT_DOUBLE_EQ(UniformBase::log_pdf(-1.423, min, max), std::log(0.2)); + EXPECT_DOUBLE_EQ(UniformBase::log_pdf(0., min, max), std::log(0.2)); + EXPECT_DOUBLE_EQ(UniformBase::log_pdf(1.31, min, max), std::log(0.2)); + EXPECT_DOUBLE_EQ(UniformBase::log_pdf(2.41, min, max), std::log(0.2)); + EXPECT_DOUBLE_EQ(UniformBase::log_pdf(2.69999999999, min, max), std::log(0.2)); +} + +TEST_F(density_fixture, uniform_log_pdf_out_of_range) +{ + EXPECT_DOUBLE_EQ(UniformBase::log_pdf(-100, min, max), std::numeric_limits::lowest()); + EXPECT_DOUBLE_EQ(UniformBase::log_pdf(-3.41, min, max), std::numeric_limits::lowest()); + EXPECT_DOUBLE_EQ(UniformBase::log_pdf(-2.3, min, max), std::numeric_limits::lowest()); + EXPECT_DOUBLE_EQ(UniformBase::log_pdf(2.7, min, max), std::numeric_limits::lowest()); + EXPECT_DOUBLE_EQ(UniformBase::log_pdf(3.5, min, max), std::numeric_limits::lowest()); + EXPECT_DOUBLE_EQ(UniformBase::log_pdf(3214, min, max), std::numeric_limits::lowest()); +} + +TEST_F(density_fixture, normal_pdf) +{ + EXPECT_NEAR(NormalBase::pdf(-10.231, mean, stddev), 1.726752595588348216742E-15, tol); + EXPECT_NEAR(NormalBase::pdf(-5.31, mean, stddev), 2.774166877919518907166E-5, tol); + EXPECT_DOUBLE_EQ(NormalBase::pdf(-2.3141231, mean, stddev), 0.04063645713784323551341); + EXPECT_DOUBLE_EQ(NormalBase::pdf(0., mean, stddev), 0.2988151821496727914542); + EXPECT_DOUBLE_EQ(NormalBase::pdf(1.31, mean, stddev), 0.2269313951019926611687); + EXPECT_DOUBLE_EQ(NormalBase::pdf(3.21, mean, stddev), 0.02505560241243631472997); + EXPECT_NEAR(NormalBase::pdf(5.24551, mean, stddev), 2.20984513448306056291E-4, tol); + EXPECT_NEAR(NormalBase::pdf(10.5699, mean, stddev), 8.61135160183067521907E-15, tol); +} + +TEST_F(density_fixture, normal_log_pdf) +{ + EXPECT_DOUBLE_EQ(NormalBase::log_pdf(-10.231, mean, stddev), std::log(1.726752595588348216742E-15)); + EXPECT_DOUBLE_EQ(NormalBase::log_pdf(-5.31, mean, stddev), std::log(2.774166877919518907166E-5)); + EXPECT_DOUBLE_EQ(NormalBase::log_pdf(-2.3141231, mean, stddev), std::log(0.04063645713784323551341)); + EXPECT_DOUBLE_EQ(NormalBase::log_pdf(0., mean, stddev), std::log(0.2988151821496727914542)); + EXPECT_DOUBLE_EQ(NormalBase::log_pdf(1.31, mean, stddev), std::log(0.2269313951019926611687)); + EXPECT_DOUBLE_EQ(NormalBase::log_pdf(3.21, mean, stddev), std::log(0.02505560241243631472997)); + EXPECT_DOUBLE_EQ(NormalBase::log_pdf(5.24551, mean, stddev), std::log(2.20984513448306056291E-4)); + EXPECT_DOUBLE_EQ(NormalBase::log_pdf(10.5699, mean, stddev), std::log(8.61135160183067521907E-15)); +} + +/* + * Discrete distributions + */ + +TEST_F(density_fixture, bernoulli_pdf_in_range) +{ + EXPECT_DOUBLE_EQ(BernoulliBase::pdf(0, p), 1-p); + EXPECT_DOUBLE_EQ(BernoulliBase::pdf(1, p), p); +} + +TEST_F(density_fixture, bernoulli_pdf_out_of_range) +{ + EXPECT_DOUBLE_EQ(BernoulliBase::pdf(-100, p), 0.); + EXPECT_DOUBLE_EQ(BernoulliBase::pdf(-3.41, p), 0.); + EXPECT_DOUBLE_EQ(BernoulliBase::pdf(-0.00000001, p), 0.); + EXPECT_DOUBLE_EQ(BernoulliBase::pdf(0.00000001, p), 0.); + EXPECT_DOUBLE_EQ(BernoulliBase::pdf(0.99999999, p), 0.); + EXPECT_DOUBLE_EQ(BernoulliBase::pdf(1.00000001, p), 0.); + EXPECT_DOUBLE_EQ(BernoulliBase::pdf(3.1423, p), 0.); + EXPECT_DOUBLE_EQ(BernoulliBase::pdf(5.613, p), 0.); + EXPECT_DOUBLE_EQ(BernoulliBase::pdf(100., p), 0.); +} + +TEST_F(density_fixture, bernoulli_pdf_always_tail) +{ + EXPECT_DOUBLE_EQ(BernoulliBase::pdf(0, 0.), 1.); + EXPECT_DOUBLE_EQ(BernoulliBase::pdf(1, 0.), 0.); +} + +TEST_F(density_fixture, bernoulli_pdf_always_head) +{ + EXPECT_DOUBLE_EQ(BernoulliBase::pdf(0, 1.), 0.); + EXPECT_DOUBLE_EQ(BernoulliBase::pdf(1, 1.), 1.); +} + +TEST_F(density_fixture, bernoulli_log_pdf_in_range) +{ + EXPECT_DOUBLE_EQ(BernoulliBase::log_pdf(0, p), std::log(1-p)); + EXPECT_DOUBLE_EQ(BernoulliBase::log_pdf(1, p), std::log(p)); +} + +TEST_F(density_fixture, bernoulli_log_pdf_out_of_range) +{ + EXPECT_DOUBLE_EQ(BernoulliBase::log_pdf(-100, p), std::numeric_limits::lowest()); + EXPECT_DOUBLE_EQ(BernoulliBase::log_pdf(-3.41, p), std::numeric_limits::lowest()); + EXPECT_DOUBLE_EQ(BernoulliBase::log_pdf(-0.00000001, p), std::numeric_limits::lowest()); + EXPECT_DOUBLE_EQ(BernoulliBase::log_pdf(0.00000001, p), std::numeric_limits::lowest()); + EXPECT_DOUBLE_EQ(BernoulliBase::log_pdf(0.99999999, p), std::numeric_limits::lowest()); + EXPECT_DOUBLE_EQ(BernoulliBase::log_pdf(1.00000001, p), std::numeric_limits::lowest()); + EXPECT_DOUBLE_EQ(BernoulliBase::log_pdf(3.1423, p), std::numeric_limits::lowest()); + EXPECT_DOUBLE_EQ(BernoulliBase::log_pdf(5.613, p), std::numeric_limits::lowest()); + EXPECT_DOUBLE_EQ(BernoulliBase::log_pdf(100., p), std::numeric_limits::lowest()); +} + +TEST_F(density_fixture, bernoulli_log_pdf_always_tail) +{ + EXPECT_DOUBLE_EQ(BernoulliBase::log_pdf(0, 0.), 0.); + EXPECT_DOUBLE_EQ(BernoulliBase::log_pdf(1, 0.), std::numeric_limits::lowest()); +} + +TEST_F(density_fixture, bernoulli_log_pdf_always_head) +{ + EXPECT_DOUBLE_EQ(BernoulliBase::log_pdf(0, 1.), std::numeric_limits::lowest()); + EXPECT_DOUBLE_EQ(BernoulliBase::log_pdf(1, 1.), 0.); +} + +} // namespace expr +} // namespace ppl diff --git a/test/expression/distribution/normal_unittest.cpp b/test/expression/distribution/normal_unittest.cpp new file mode 100644 index 00000000..27b12291 --- /dev/null +++ b/test/expression/distribution/normal_unittest.cpp @@ -0,0 +1,68 @@ +#include "gtest/gtest.h" +#include +#include +#include +#include +#include + +namespace ppl { +namespace expr { + +struct normal_fixture : ::testing::Test { +protected: + using value_t = typename MockVarExpr::value_t; + static constexpr size_t sample_size = 1000; + double mean = 0.1; + double stddev = 0.8; + MockVarExpr x{mean}; + MockVarExpr y{stddev}; + using norm_t = Normal; + norm_t norm = {x, y}; + std::array sample = {0.}; +}; + +TEST_F(normal_fixture, ctor) +{ + static_assert(util::is_dist_expr_v); +} + +TEST_F(normal_fixture, normal_check_params) { + EXPECT_DOUBLE_EQ(norm.mean(), static_cast(x)); + EXPECT_DOUBLE_EQ(norm.stddev(), static_cast(y)); +} + +TEST_F(normal_fixture, normal_pdf_delegate) { + EXPECT_DOUBLE_EQ(norm.pdf(-10.664), NormalBase::pdf(-10.664, mean, stddev)); + EXPECT_DOUBLE_EQ(norm.pdf(-7.324), NormalBase::pdf(-7.324, mean, stddev)); + EXPECT_DOUBLE_EQ(norm.pdf(-3.241), NormalBase::pdf(-3.241, mean, stddev)); + EXPECT_DOUBLE_EQ(norm.pdf(-0.359288), NormalBase::pdf(-0.359288, mean, stddev)); + EXPECT_DOUBLE_EQ(norm.pdf(0.12314), NormalBase::pdf(0.12314, mean, stddev)); + EXPECT_DOUBLE_EQ(norm.pdf(3.145), NormalBase::pdf(3.145, mean, stddev)); + EXPECT_DOUBLE_EQ(norm.pdf(6.000923), NormalBase::pdf(6.000923, mean, stddev)); + EXPECT_DOUBLE_EQ(norm.pdf(16.423), NormalBase::pdf(16.423, mean, stddev)); +} + +TEST_F(normal_fixture, normal_log_pdf_delegate) { + EXPECT_DOUBLE_EQ(norm.log_pdf(-10.664), NormalBase::log_pdf(-10.664, mean, stddev)); + EXPECT_DOUBLE_EQ(norm.log_pdf(-7.324), NormalBase::log_pdf(-7.324, mean, stddev)); + EXPECT_DOUBLE_EQ(norm.log_pdf(-3.241), NormalBase::log_pdf(-3.241, mean, stddev)); + EXPECT_DOUBLE_EQ(norm.log_pdf(-0.359288), NormalBase::log_pdf(-0.359288, mean, stddev)); + EXPECT_DOUBLE_EQ(norm.log_pdf(0.12314), NormalBase::log_pdf(0.12314, mean, stddev)); + EXPECT_DOUBLE_EQ(norm.log_pdf(3.145), NormalBase::log_pdf(3.145, mean, stddev)); + EXPECT_DOUBLE_EQ(norm.log_pdf(6.000923), NormalBase::log_pdf(6.000923, mean, stddev)); + EXPECT_DOUBLE_EQ(norm.log_pdf(16.423), NormalBase::log_pdf(16.423, mean, stddev)); +} + +TEST_F(normal_fixture, normal_sample) { + std::random_device rd{}; + std::mt19937 gen{rd()}; + + for (size_t i = 0; i < sample_size; i++) { + sample[i] = norm.sample(gen); + } + + plot_hist(sample); +} + +} // namespace expr +} // namespace ppl diff --git a/test/expression/distribution/uniform_unittest.cpp b/test/expression/distribution/uniform_unittest.cpp new file mode 100644 index 00000000..7177ffb0 --- /dev/null +++ b/test/expression/distribution/uniform_unittest.cpp @@ -0,0 +1,70 @@ +#include "gtest/gtest.h" +#include +#include +#include +#include +#include + +namespace ppl { +namespace expr { + +struct uniform_fixture : ::testing::Test { +protected: + using value_t = typename MockVarExpr::value_t; + static constexpr size_t sample_size = 1000; + double min = 0.1; + double max = 0.8; + MockVarExpr x{min}; + MockVarExpr y{max}; + using unif_t = Uniform; + unif_t unif = {x, y}; + std::array sample = {0.}; +}; + +TEST_F(uniform_fixture, ctor) +{ + static_assert(util::is_dist_expr_v); +} + +TEST_F(uniform_fixture, uniform_check_params) { + EXPECT_DOUBLE_EQ(unif.min(), static_cast(x)); + EXPECT_DOUBLE_EQ(unif.max(), static_cast(y)); +} + +TEST_F(uniform_fixture, uniform_pdf_delegate) { + EXPECT_DOUBLE_EQ(unif.pdf(-10.664), UniformBase::pdf(-10.664, min, max)); + EXPECT_DOUBLE_EQ(unif.pdf(-7.324), UniformBase::pdf(-7.324, min, max)); + EXPECT_DOUBLE_EQ(unif.pdf(-3.241), UniformBase::pdf(-3.241, min, max)); + EXPECT_DOUBLE_EQ(unif.pdf(-0.359288), UniformBase::pdf(-0.359288, min, max)); + EXPECT_DOUBLE_EQ(unif.pdf(0.12314), UniformBase::pdf(0.12314, min, max)); + EXPECT_DOUBLE_EQ(unif.pdf(3.145), UniformBase::pdf(3.145, min, max)); + EXPECT_DOUBLE_EQ(unif.pdf(6.000923), UniformBase::pdf(6.000923, min, max)); + EXPECT_DOUBLE_EQ(unif.pdf(16.423), UniformBase::pdf(16.423, min, max)); +} + +TEST_F(uniform_fixture, uniform_log_pdf_delegate) { + EXPECT_DOUBLE_EQ(unif.log_pdf(-10.664), UniformBase::log_pdf(-10.664, min, max)); + EXPECT_DOUBLE_EQ(unif.log_pdf(-7.324), UniformBase::log_pdf(-7.324, min, max)); + EXPECT_DOUBLE_EQ(unif.log_pdf(-3.241), UniformBase::log_pdf(-3.241, min, max)); + EXPECT_DOUBLE_EQ(unif.log_pdf(-0.359288), UniformBase::log_pdf(-0.359288, min, max)); + EXPECT_DOUBLE_EQ(unif.log_pdf(0.12314), UniformBase::log_pdf(0.12314, min, max)); + EXPECT_DOUBLE_EQ(unif.log_pdf(3.145), UniformBase::log_pdf(3.145, min, max)); + EXPECT_DOUBLE_EQ(unif.log_pdf(6.000923), UniformBase::log_pdf(6.000923, min, max)); + EXPECT_DOUBLE_EQ(unif.log_pdf(16.423), UniformBase::log_pdf(16.423, min, max)); +} + +TEST_F(uniform_fixture, uniform_sample) { + std::random_device rd{}; + std::mt19937 gen{rd()}; + + for (size_t i = 0; i < sample_size; i++) { + sample[i] = unif.sample(gen); + EXPECT_GT(sample[i], min); + EXPECT_LT(sample[i], max); + } + + plot_hist(sample, 0.05, min, max); +} + +} // namespace expr +} // namespace ppl diff --git a/test/expression/model_unittest.cpp b/test/expression/model/model_unittest.cpp similarity index 53% rename from test/expression/model_unittest.cpp rename to test/expression/model/model_unittest.cpp index 4cc39eec..b8ab9625 100644 --- a/test/expression/model_unittest.cpp +++ b/test/expression/model/model_unittest.cpp @@ -1,7 +1,8 @@ +#include "gtest/gtest.h" #include #include -#include "gtest/gtest.h" -#include +#include +#include namespace ppl { namespace expr { @@ -10,55 +11,6 @@ namespace expr { // Model with one RV TESTS ////////////////////////////////////////////////////// -/* - * Mock state class for testing purposes. - */ -enum class MockState { - data, - parameter -}; - -/* - * Mock var object for testing purposes. - * Must meet some of the requirements of actual var types. - */ -struct MockVar -{ - using value_t = double; - using pointer_t = double*; - using state_t = MockState; - - void set_value(double val) { value_ = val; } - value_t get_value() const { return value_; } - - void set_state(state_t state) { state_ = state; } - state_t get_state() const { return state_; } - - operator value_t() { return value_; } - -private: - double value_; - state_t state_ = state_t::parameter; -}; - -/* - * Mock distribution object for testing purposes. - * Must meet some of the requirements of actual distribution types. - */ -struct MockDist -{ - using value_t = double; - using dist_value_t = double; - - // Mock pdf - identity function. - double pdf(double x) const - { return x; } - - // Mock log_pdf - log(pdf(x)). - double log_pdf(double x) const - { return std::log(pdf(x)); } -}; - /* * Fixture for testing one var with distribution. */ @@ -66,84 +18,49 @@ struct var_dist_fixture : ::testing::Test { protected: MockVar x; - using model_t = EqNode; - model_t model = {x, MockDist()}; + using model_t = EqNode; + model_t model = {x, MockDistExpr()}; double val; - void reconfigure(double val) - { - x.set_value(val); - } + void reconfigure() + { x.set_value(val); } }; +TEST_F(var_dist_fixture, ctor) +{ + static_assert(util::is_model_expr_v); +} + TEST_F(var_dist_fixture, pdf_valid) { - // MockDist pdf is identity function + // MockDistExpr pdf is identity function // so we may simply compare model.pdf() with val. val = 0.000001; - reconfigure(val); + reconfigure(); EXPECT_EQ(model.pdf(), val); val = 0.5; - reconfigure(val); + reconfigure(); EXPECT_EQ(model.pdf(), val); val = 0.999999; - reconfigure(val); - EXPECT_EQ(model.pdf(), val); -} - -TEST_F(var_dist_fixture, pdf_invalid) -{ - val = 0.000004123; - reconfigure(val); - EXPECT_EQ(model.pdf(), val); - - val = 0.55555; - reconfigure(val); - EXPECT_EQ(model.pdf(), val); - - val = 5.234424231; - reconfigure(val); - EXPECT_EQ(model.pdf(), val); - - val = 69; - reconfigure(val); + reconfigure(); EXPECT_EQ(model.pdf(), val); } TEST_F(var_dist_fixture, log_pdf_valid) { val = 0.000001; - reconfigure(val); + reconfigure(); EXPECT_EQ(model.log_pdf(), std::log(val)); val = 0.5; - reconfigure(val); + reconfigure(); EXPECT_EQ(model.log_pdf(), std::log(val)); val = 0.999999; - reconfigure(val); - EXPECT_EQ(model.log_pdf(), std::log(val)); -} - -TEST_F(var_dist_fixture, log_pdf_invalid) -{ - val = 0.000004123; - reconfigure(val); - EXPECT_EQ(model.log_pdf(), std::log(val)); - - val = 0.55555; - reconfigure(val); - EXPECT_EQ(model.log_pdf(), std::log(val)); - - val = 5.234424231; - reconfigure(val); - EXPECT_EQ(model.log_pdf(), std::log(val)); - - val = 69; - reconfigure(val); + reconfigure(); EXPECT_EQ(model.log_pdf(), std::log(val)); } @@ -159,13 +76,12 @@ struct many_var_dist_fixture : ::testing::Test protected: MockVar x, y, z, w; double xv, yv, zv, wv; - using eq_t = EqNode; + using eq_t = EqNode; using model_two_t = GlueNode; model_two_t model_two = { - {x, MockDist()}, - {y, MockDist()} - + {x, MockDistExpr()}, + {y, MockDistExpr()} }; using model_four_t = @@ -176,17 +92,23 @@ struct many_var_dist_fixture : ::testing::Test >; model_four_t model_four = { - {x, MockDist()}, + {x, MockDistExpr()}, { - {y, MockDist()}, + {y, MockDistExpr()}, { - {z, MockDist()}, - {w, MockDist()} + {z, MockDistExpr()}, + {w, MockDistExpr()} } } }; }; +TEST_F(many_var_dist_fixture, ctor) +{ + static_assert(util::is_model_expr_v); + static_assert(util::is_model_expr_v); +} + TEST_F(many_var_dist_fixture, two_vars_pdf) { xv = 0.2; yv = 1.8; @@ -218,7 +140,7 @@ TEST_F(many_var_dist_fixture, four_vars_traverse_count_params) z.set_state(MockState::data); model_four.traverse([&](auto& model) { using var_t = std::decay_t; - using state_t = typename var_traits::state_t; + using state_t = typename util::var_traits::state_t; count += (model.get_variable().get_state() == state_t::parameter); }); EXPECT_EQ(count, 3); diff --git a/test/expression/variable/constant_unittest.cpp b/test/expression/variable/constant_unittest.cpp new file mode 100644 index 00000000..6efdc71e --- /dev/null +++ b/test/expression/variable/constant_unittest.cpp @@ -0,0 +1,30 @@ +#include "gtest/gtest.h" +#include +#include +#include + +namespace ppl { +namespace expr { + +struct constant_fixture : ::testing::Test +{ +protected: + using value_t = double; + value_t c = 0.3; + Constant x{c}; +}; + +TEST_F(constant_fixture, ctor) +{ + static_assert(util::is_var_expr_v>); +} + +TEST_F(constant_fixture, convertible_value) +{ + EXPECT_EQ(static_cast(x), 0.3); + c = 3.41; + EXPECT_EQ(static_cast(x), 0.3); +} + +} // namespace expr +} // namespace ppl diff --git a/test/expression/variable/variable_viewer_unittest.cpp b/test/expression/variable/variable_viewer_unittest.cpp new file mode 100644 index 00000000..1324efed --- /dev/null +++ b/test/expression/variable/variable_viewer_unittest.cpp @@ -0,0 +1,32 @@ +#include "gtest/gtest.h" +#include +#include + +namespace ppl { +namespace expr { + +struct variable_viewer_fixture : ::testing::Test +{ +protected: + using value_t = typename MockVar::value_t; + MockVar var; + VariableViewer x = var; +}; + +TEST_F(variable_viewer_fixture, ctor) +{ + static_assert(util::is_var_expr_v>); +} + +TEST_F(variable_viewer_fixture, convertible_value) +{ + var.set_value(1.); + EXPECT_EQ(static_cast(x), 1.); + + // Tests if viewer correctly reflects any changes that happened in var. + var.set_value(-3.14); + EXPECT_EQ(static_cast(x), -3.14); +} + +} // namespace expr +} // namespace ppl diff --git a/test/testutil/mock_types.hpp b/test/testutil/mock_types.hpp new file mode 100644 index 00000000..fe4cbee1 --- /dev/null +++ b/test/testutil/mock_types.hpp @@ -0,0 +1,118 @@ +#pragma once +#include + +namespace ppl { + +/* + * Mock state class for testing purposes. + */ +enum class MockState { + data, + parameter +}; + +/* + * Mock Variable class that should meet the requirements + * of is_var_v. + */ +struct MockVar +{ + using value_t = double; + using pointer_t = double*; + using const_pointer_t = const double*; + using state_t = MockState; + + operator value_t() const {return x_;} + void set_value(value_t x) {x_ = x;} + value_t get_value() const {return x_;} + + void set_storage(pointer_t ptr) {ptr_ = ptr;} + + void set_state(state_t state) {state_ = state;} + state_t get_state() const {return state_;} + +private: + value_t x_ = 0.; + pointer_t ptr_ = nullptr; + state_t state_ = state_t::parameter; +}; + +/* + * Mock variable classes that fulfill + * var_traits requirements, but do not fit the rest. + */ +struct MockVar_no_convertible +{ + using value_t = double; + using pointer_t = double*; + using state_t = void; +}; + +/* + * Mock Variable Expression class that should meet the requirements + * of is_var_expr_v. + */ +struct MockVarExpr +{ + using value_t = double; + operator value_t() const {return x_;} + + /* not part of API */ + MockVarExpr(value_t x) + : x_{x} + {} + + void set_value(value_t x) {x_ = x;} +private: + value_t x_ = 0.; +}; + +/* + * Mock variable expression classes that fulfill + * var_expr_traits requirements, but do not fit the rest. + */ +struct MockVarExpr_no_convertible +{ + using value_t = double; +}; + +/* + * Mock distribution expression class that should meet the requirements + * of is_dist_expr_v. + */ +struct MockDistExpr +{ + using value_t = double; + using dist_value_t = double; + + dist_value_t pdf(value_t x) const + { return x; } + + dist_value_t log_pdf(value_t x) const + { return std::log(x); } +}; + +/* + * Mock distribution expression classes that fulfill + * dist_expr_traits requirements, but do not fit the rest. + */ +struct MockDistExpr_no_pdf : public MockDistExpr +{ +private: + using MockDistExpr::pdf; +}; + +struct MockDistExpr_no_log_pdf : public MockDistExpr +{ +private: + using MockDistExpr::log_pdf; +}; + +/* + * TODO: + * Mock model expression clases that should meet the + * requirements of is_model_expr_v. + * Additionally, MockEqNode should satisfy is_eq_node_expr_v. + */ + +} // namespace ppl diff --git a/test/testutil/sample_tools.hpp b/test/testutil/sample_tools.hpp new file mode 100644 index 00000000..68537d5c --- /dev/null +++ b/test/testutil/sample_tools.hpp @@ -0,0 +1,53 @@ +#pragma once +#include "gtest/gtest.h" +#include +#include +#include +#include +#include +#include + +namespace ppl { + +// Plotting utility to visualize histogram of samples. +template +inline void plot_hist(const ArrayType& arr, + double step_size = .5, + double min = std::numeric_limits::lowest(), + double max = std::numeric_limits::max() + ) +{ + constexpr size_t nstars = 100; // maximum number of stars to distribute + + min = (min == std::numeric_limits::lowest()) ? + *std::min_element(arr.begin(), arr.end()) : + min; + max = (max == std::numeric_limits::max()) ? + *std::max_element(arr.begin(), arr.end()) : + max; + const int64_t nearest_min = std::floor(min); + const int64_t nearest_max = std::floor(max) + 1; + const uint64_t range = nearest_max - nearest_min; + const uint64_t n_hist = std::floor(range/step_size); + + // keeps count for each histogram bar + std::vector counter(n_hist, 0); + + for (auto x : arr) { + if (nearest_min <= x && x <= nearest_max) { + ++counter[std::floor((x - nearest_min) / step_size)]; + } + } + + if ((min == *std::min_element(arr.begin(), arr.end())) && + (max == *std::max_element(arr.begin(), arr.end()))) { + EXPECT_EQ(std::accumulate(counter.begin(), counter.end(), 0), (int) arr.size()); + } + + for (size_t i = 0; i < n_hist; ++i) { + std::cout << i << "-" << (i+1) << ": " << '\t'; + std::cout << std::string(counter[i] * nstars/arr.size(), '*') << std::endl; + } +} + +} // namespace ppl diff --git a/test/util/concept_unittest.cpp b/test/util/concept_unittest.cpp new file mode 100644 index 00000000..d2fe4899 --- /dev/null +++ b/test/util/concept_unittest.cpp @@ -0,0 +1,91 @@ +#include "gtest/gtest.h" +#include + +namespace ppl { +namespace util { + +struct MockType +{}; + +struct MockType2 +{ + using value_t = double; + using pointer_t = double*; + using state_t = void; + + void pdf() {}; + void log_pdf() {}; +}; + +struct concept_fixture : ::testing::Test +{ +protected: +}; + +TEST_F(concept_fixture, has_type_value_t_v_false) +{ + static_assert(!has_type_value_t_v); + static_assert(!has_type_value_t_v); + static_assert(!has_type_value_t_v); + static_assert(!has_type_value_t_v); +} + +TEST_F(concept_fixture, has_type_value_t_v_true) +{ + static_assert(has_type_value_t_v); +} + +TEST_F(concept_fixture, has_type_pointer_t_v_false) +{ + static_assert(!has_type_pointer_t_v); + static_assert(!has_type_pointer_t_v); + static_assert(!has_type_pointer_t_v); + static_assert(!has_type_pointer_t_v); +} + +TEST_F(concept_fixture, has_type_pointer_t_v_true) +{ + static_assert(has_type_pointer_t_v); +} + +TEST_F(concept_fixture, has_type_state_t_v_false) +{ + static_assert(!has_type_state_t_v); + static_assert(!has_type_state_t_v); + static_assert(!has_type_state_t_v); + static_assert(!has_type_state_t_v); +} + +TEST_F(concept_fixture, has_type_state_t_v_true) +{ + static_assert(has_type_state_t_v); +} + +TEST_F(concept_fixture, has_func_pdf_v_false) +{ + static_assert(!has_func_pdf_v); + static_assert(!has_func_pdf_v); + static_assert(!has_func_pdf_v); + static_assert(!has_func_pdf_v); +} + +TEST_F(concept_fixture, has_func_pdf_v_true) +{ + static_assert(has_func_pdf_v); +} + +TEST_F(concept_fixture, has_func_log_pdf_v_false) +{ + static_assert(!has_func_log_pdf_v); + static_assert(!has_func_log_pdf_v); + static_assert(!has_func_log_pdf_v); + static_assert(!has_func_log_pdf_v); +} + +TEST_F(concept_fixture, has_func_log_pdf_v_true) +{ + static_assert(has_func_log_pdf_v); +} + +} // namespace util +} // namespace ppl diff --git a/test/util/dist_expr_traits_unittest.cpp b/test/util/dist_expr_traits_unittest.cpp new file mode 100644 index 00000000..3858405f --- /dev/null +++ b/test/util/dist_expr_traits_unittest.cpp @@ -0,0 +1,25 @@ +#include "gtest/gtest.h" +#include +#include + +namespace ppl { +namespace util { + +struct dist_expr_traits_fixture : ::testing::Test +{ +protected: +}; + +TEST_F(dist_expr_traits_fixture, is_dist_expr_v_true) +{ + static_assert(is_dist_expr_v); +} + +TEST_F(dist_expr_traits_fixture, is_dist_expr_v_false) +{ + static_assert(!is_dist_expr_v); + static_assert(!is_dist_expr_v); +} + +} // namespace util +} // namespace ppl diff --git a/test/util/var_expr_traits_unittest.cpp b/test/util/var_expr_traits_unittest.cpp new file mode 100644 index 00000000..9ca6d875 --- /dev/null +++ b/test/util/var_expr_traits_unittest.cpp @@ -0,0 +1,24 @@ +#include "gtest/gtest.h" +#include +#include + +namespace ppl { +namespace util { + +struct var_expr_traits_fixture : ::testing::Test +{ +protected: +}; + +TEST_F(var_expr_traits_fixture, is_var_expr_v_true) +{ + static_assert(is_var_expr_v); +} + +TEST_F(var_expr_traits_fixture, is_var_expr_v_false) +{ + static_assert(!is_var_expr_v); +} + +} // namespace util +} // namespace ppl diff --git a/test/util/var_traits_unittest.cpp b/test/util/var_traits_unittest.cpp new file mode 100644 index 00000000..30823073 --- /dev/null +++ b/test/util/var_traits_unittest.cpp @@ -0,0 +1,24 @@ +#include "gtest/gtest.h" +#include +#include + +namespace ppl { +namespace util { + +struct var_traits_fixture : ::testing::Test +{ +protected: +}; + +TEST_F(var_traits_fixture, is_var_v_true) +{ + static_assert(is_var_v); +} + +TEST_F(var_traits_fixture, is_var_v_false) +{ + static_assert(!is_var_v); +} + +} // namespace util +} // namespace ppl From 7af98347602166ac4caaf8a1f110bd79c7f221d2 Mon Sep 17 00:00:00 2001 From: James Yang Date: Mon, 20 Apr 2020 17:56:05 -0400 Subject: [PATCH 33/35] Add TODOs and extend tolerance of sample mean --- include/autoppl/expr_builder.hpp | 1 + test/algorithm/mh_unittest.cpp | 18 +++--------------- 2 files changed, 4 insertions(+), 15 deletions(-) diff --git a/include/autoppl/expr_builder.hpp b/include/autoppl/expr_builder.hpp index 53f50449..43365307 100644 --- a/include/autoppl/expr_builder.hpp +++ b/include/autoppl/expr_builder.hpp @@ -117,6 +117,7 @@ inline constexpr auto normal(MeanType&& mean_expr, * Builds a Bernoulli expression only when the parameter * is a valid discrete distribution parameter type. * See var_expr.hpp for more information. + * TODO: generalize as done with uniform and normal */ template inline constexpr auto bernoulli(const ProbType& p_expr) diff --git a/test/algorithm/mh_unittest.cpp b/test/algorithm/mh_unittest.cpp index b88400ac..c1fd9c20 100644 --- a/test/algorithm/mh_unittest.cpp +++ b/test/algorithm/mh_unittest.cpp @@ -7,18 +7,6 @@ namespace ppl { -TEST(mh_unittest, plot_hist_sanity) -{ - static constexpr size_t sample_size = 20000; - std::array storage = {0.}; - std::normal_distribution normal_sampler(0., 1.); - std::mt19937 gen; - for (size_t i = 0; i < sample_size; ++i) { - storage[i] = normal_sampler(gen); - } - plot_hist(storage); -} - /* * Fixture for Metropolis-Hastings */ @@ -48,7 +36,7 @@ TEST_F(mh_fixture, sample_std_normal) auto model = (theta |= normal(0., 1.)); mh_posterior(model, sample_size); plot_hist(storage); - EXPECT_NEAR(sample_average(1000), 0., 0.01); + EXPECT_NEAR(sample_average(1000), 0., 0.1); } TEST_F(mh_fixture, sample_uniform) @@ -56,7 +44,7 @@ TEST_F(mh_fixture, sample_uniform) auto model = (theta |= uniform(0., 1.)); mh_posterior(model, sample_size); plot_hist(storage, 0.05, 0., 1.); - EXPECT_NEAR(sample_average(1000), 0.5, 0.01); + EXPECT_NEAR(sample_average(1000), 0.5, 0.1); } TEST_F(mh_fixture, sample_unif_normal_posterior) @@ -68,7 +56,7 @@ TEST_F(mh_fixture, sample_unif_normal_posterior) ); mh_posterior(model, sample_size); plot_hist(storage); - EXPECT_NEAR(sample_average(1000), 3.0, 0.01); + EXPECT_NEAR(sample_average(1000), 3.0, 0.1); } } // namespace ppl From 1cc7b9f6ea44309f336123013694c12b3fc0c83a Mon Sep 17 00:00:00 2001 From: James Yang Date: Mon, 20 Apr 2020 19:20:02 -0400 Subject: [PATCH 34/35] MH default seed to current time and tests use fixed seed of 0 --- include/autoppl/algorithm/mh.hpp | 3 ++- test/algorithm/mh_unittest.cpp | 15 ++++++++------- 2 files changed, 10 insertions(+), 8 deletions(-) diff --git a/include/autoppl/algorithm/mh.hpp b/include/autoppl/algorithm/mh.hpp index 1a83d3f8..d9872524 100644 --- a/include/autoppl/algorithm/mh.hpp +++ b/include/autoppl/algorithm/mh.hpp @@ -1,4 +1,5 @@ #pragma once +#include #include #include #include @@ -36,7 +37,7 @@ template inline void mh_posterior(ModelType& model, double n_sample, double stddev = 1.0, - double seed = 0 + double seed = std::time(0) ) { using data_t = details::MHData; diff --git a/test/algorithm/mh_unittest.cpp b/test/algorithm/mh_unittest.cpp index c1fd9c20..0922112b 100644 --- a/test/algorithm/mh_unittest.cpp +++ b/test/algorithm/mh_unittest.cpp @@ -16,12 +16,13 @@ struct mh_fixture : ::testing::Test static constexpr size_t sample_size = 20000; std::array storage = {0.}; Variable theta, x; + size_t burn = 1000; mh_fixture() : theta{storage.data()} {} - double sample_average(size_t burn) + double sample_average() { double sum = std::accumulate( std::next(storage.begin(), burn), @@ -34,17 +35,17 @@ struct mh_fixture : ::testing::Test TEST_F(mh_fixture, sample_std_normal) { auto model = (theta |= normal(0., 1.)); - mh_posterior(model, sample_size); + mh_posterior(model, sample_size, 1.0, 0.); plot_hist(storage); - EXPECT_NEAR(sample_average(1000), 0., 0.1); + EXPECT_NEAR(sample_average(), 0., 0.1); } TEST_F(mh_fixture, sample_uniform) { auto model = (theta |= uniform(0., 1.)); - mh_posterior(model, sample_size); + mh_posterior(model, sample_size, 1.0, 0.); plot_hist(storage, 0.05, 0., 1.); - EXPECT_NEAR(sample_average(1000), 0.5, 0.1); + EXPECT_NEAR(sample_average(), 0.5, 0.1); } TEST_F(mh_fixture, sample_unif_normal_posterior) @@ -54,9 +55,9 @@ TEST_F(mh_fixture, sample_unif_normal_posterior) theta |= uniform(-20., 20.), x |= normal(theta, 1.) ); - mh_posterior(model, sample_size); + mh_posterior(model, sample_size, 1.0, 0.); plot_hist(storage); - EXPECT_NEAR(sample_average(1000), 3.0, 0.1); + EXPECT_NEAR(sample_average(), 3.0, 0.1); } } // namespace ppl From d2a1f4d0b9b5daab9fc58b0f5c7ffec97f370ace Mon Sep 17 00:00:00 2001 From: James Yang Date: Mon, 20 Apr 2020 19:29:35 -0400 Subject: [PATCH 35/35] MH default seed using chrono to get milliseconds since epoch --- include/autoppl/algorithm/mh.hpp | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/include/autoppl/algorithm/mh.hpp b/include/autoppl/algorithm/mh.hpp index d9872524..f5f17fa1 100644 --- a/include/autoppl/algorithm/mh.hpp +++ b/include/autoppl/algorithm/mh.hpp @@ -1,5 +1,5 @@ #pragma once -#include +#include #include #include #include @@ -37,8 +37,11 @@ template inline void mh_posterior(ModelType& model, double n_sample, double stddev = 1.0, - double seed = std::time(0) - ) + double seed = std::chrono::duration_cast< + std::chrono::milliseconds>( + std::chrono::system_clock::now().time_since_epoch() + ).count() + ) { using data_t = details::MHData;