diff --git a/sample/sample.cpp b/sample/sample.cpp index b505e12..53b629d 100644 --- a/sample/sample.cpp +++ b/sample/sample.cpp @@ -9,6 +9,7 @@ #include #include #include +#include #include "Models.h" #include "modelRegister.h" #include "config_platform.h" @@ -160,82 +161,119 @@ int main(int argc, char** argv) states[feature] = std::vector(maxes[feature]); } states[className] = std::vector(maxes[className]); - auto clf = platform::Models::instance()->create(model_name); + // Output the states + std::cout << std::string(80, '-') << std::endl; + std::cout << "States" << std::endl; + for (auto feature : features) { + std::cout << feature << ": " << states[feature].size() << std::endl; + } + std::cout << std::string(80, '-') << std::endl; + //auto clf = platform::Models::instance()->create("SPODE"); + auto clf = bayesnet::SPODE(2); + bayesnet::Smoothing_t smoothing = bayesnet::Smoothing_t::ORIGINAL; - clf->fit(Xd, y, features, className, states, smoothing); + clf.fit(Xd, y, features, className, states, smoothing); if (dump_cpt) { std::cout << "--- CPT Tables ---" << std::endl; - clf->dump_cpt(); + std::cout << clf.dump_cpt(); } - auto lines = clf->show(); + std::cout << "--- Datos predicción ---" << std::endl; + std::cout << "Orden de variables: " << std::endl; + for (auto feature : features) { + std::cout << feature << ", "; + } + std::cout << std::endl; + std::cout << "X[0]: "; + for (int i = 0; i < Xd.size(); ++i) { + std::cout << Xd[i][0] << ", "; + } + std::cout << std::endl; + std::cout << std::string(80, '-') << std::endl; + + auto lines = clf.show(); for (auto line : lines) { std::cout << line << std::endl; } std::cout << "--- Topological Order ---" << std::endl; - auto order = clf->topological_order(); + auto order = clf.topological_order(); for (auto name : order) { std::cout << name << ", "; } - std::cout << "end." << std::endl; - auto score = clf->score(Xd, y); - std::cout << "Score: " << score << std::endl; - auto graph = clf->graph(); - auto dot_file = model_name + "_" + file_name; - ofstream file(dot_file + ".dot"); - file << graph; - file.close(); - std::cout << "Graph saved in " << model_name << "_" << file_name << ".dot" << std::endl; - std::cout << "dot -Tpng -o " + dot_file + ".png " + dot_file + ".dot " << std::endl; - std::string stratified_string = stratified ? " Stratified" : ""; - std::cout << nFolds << " Folds" << stratified_string << " Cross validation" << std::endl; - std::cout << "==========================================" << std::endl; - torch::Tensor Xt = torch::zeros({ static_cast(Xd.size()), static_cast(Xd[0].size()) }, torch::kInt32); - torch::Tensor yt = torch::tensor(y, torch::kInt32); - for (int i = 0; i < features.size(); ++i) { - Xt.index_put_({ i, "..." }, torch::tensor(Xd[i], torch::kInt32)); - } - float total_score = 0, total_score_train = 0, score_train, score_test; - folding::Fold* fold; - double nodes = 0.0; - if (stratified) - fold = new folding::StratifiedKFold(nFolds, y, seed); - else - fold = new folding::KFold(nFolds, y.size(), seed); - for (auto i = 0; i < nFolds; ++i) { - auto [train, test] = fold->getFold(i); - std::cout << "Fold: " << i + 1 << std::endl; - if (tensors) { - auto ttrain = torch::tensor(train, torch::kInt64); - auto ttest = torch::tensor(test, torch::kInt64); - torch::Tensor Xtraint = torch::index_select(Xt, 1, ttrain); - torch::Tensor ytraint = yt.index({ ttrain }); - torch::Tensor Xtestt = torch::index_select(Xt, 1, ttest); - torch::Tensor ytestt = yt.index({ ttest }); - clf->fit(Xtraint, ytraint, features, className, states, smoothing); - auto temp = clf->predict(Xtraint); - score_train = clf->score(Xtraint, ytraint); - score_test = clf->score(Xtestt, ytestt); - } else { - auto [Xtrain, ytrain] = extract_indices(train, Xd, y); - auto [Xtest, ytest] = extract_indices(test, Xd, y); - clf->fit(Xtrain, ytrain, features, className, states, smoothing); - std::cout << "Nodes: " << clf->getNumberOfNodes() << std::endl; - nodes += clf->getNumberOfNodes(); - score_train = clf->score(Xtrain, ytrain); - score_test = clf->score(Xtest, ytest); + auto predict_proba = clf.predict_proba(Xd); + std::cout << "Instances predict_proba: "; + for (int i = 0; i < predict_proba.size(); i++) { + std::cout << "Instance " << i << ": "; + for (int j = 0; j < 4; ++j) { + std::cout << Xd[j][i] << ", "; } - if (dump_cpt) { - std::cout << "--- CPT Tables ---" << std::endl; - std::cout << clf->dump_cpt(); + std::cout << ": "; + for (auto score : predict_proba[i]) { + std::cout << score << ", "; } - total_score_train += score_train; - total_score += score_test; - std::cout << "Score Train: " << score_train << std::endl; - std::cout << "Score Test : " << score_test << std::endl; - std::cout << "-------------------------------------------------------------------------------" << std::endl; + std::cout << std::endl; } - std::cout << "Nodes: " << nodes / nFolds << std::endl; - std::cout << "**********************************************************************************" << std::endl; - std::cout << "Average Score Train: " << total_score_train / nFolds << std::endl; - std::cout << "Average Score Test : " << total_score / nFolds << std::endl;return 0; + // std::cout << std::endl; + // std::cout << "end." << std::endl; + // auto score = clf->score(Xd, y); + // std::cout << "Score: " << score << std::endl; + // auto graph = clf->graph(); + // auto dot_file = model_name + "_" + file_name; + // ofstream file(dot_file + ".dot"); + // file << graph; + // file.close(); + // std::cout << "Graph saved in " << model_name << "_" << file_name << ".dot" << std::endl; + // std::cout << "dot -Tpng -o " + dot_file + ".png " + dot_file + ".dot " << std::endl; + // std::string stratified_string = stratified ? " Stratified" : ""; + // std::cout << nFolds << " Folds" << stratified_string << " Cross validation" << std::endl; + // std::cout << "==========================================" << std::endl; + // torch::Tensor Xt = torch::zeros({ static_cast(Xd.size()), static_cast(Xd[0].size()) }, torch::kInt32); + // torch::Tensor yt = torch::tensor(y, torch::kInt32); + // for (int i = 0; i < features.size(); ++i) { + // Xt.index_put_({ i, "..." }, torch::tensor(Xd[i], torch::kInt32)); + // } + // float total_score = 0, total_score_train = 0, score_train, score_test; + // folding::Fold* fold; + // double nodes = 0.0; + // if (stratified) + // fold = new folding::StratifiedKFold(nFolds, y, seed); + // else + // fold = new folding::KFold(nFolds, y.size(), seed); + // for (auto i = 0; i < nFolds; ++i) { + // auto [train, test] = fold->getFold(i); + // std::cout << "Fold: " << i + 1 << std::endl; + // if (tensors) { + // auto ttrain = torch::tensor(train, torch::kInt64); + // auto ttest = torch::tensor(test, torch::kInt64); + // torch::Tensor Xtraint = torch::index_select(Xt, 1, ttrain); + // torch::Tensor ytraint = yt.index({ ttrain }); + // torch::Tensor Xtestt = torch::index_select(Xt, 1, ttest); + // torch::Tensor ytestt = yt.index({ ttest }); + // clf->fit(Xtraint, ytraint, features, className, states, smoothing); + // auto temp = clf->predict(Xtraint); + // score_train = clf->score(Xtraint, ytraint); + // score_test = clf->score(Xtestt, ytestt); + // } else { + // auto [Xtrain, ytrain] = extract_indices(train, Xd, y); + // auto [Xtest, ytest] = extract_indices(test, Xd, y); + // clf->fit(Xtrain, ytrain, features, className, states, smoothing); + // std::cout << "Nodes: " << clf->getNumberOfNodes() << std::endl; + // nodes += clf->getNumberOfNodes(); + // score_train = clf->score(Xtrain, ytrain); + // score_test = clf->score(Xtest, ytest); + // } + // // if (dump_cpt) { + // // std::cout << "--- CPT Tables ---" << std::endl; + // // std::cout << clf->dump_cpt(); + // // } + // total_score_train += score_train; + // total_score += score_test; + // std::cout << "Score Train: " << score_train << std::endl; + // std::cout << "Score Test : " << score_test << std::endl; + // std::cout << "-------------------------------------------------------------------------------" << std::endl; + // } + + // std::cout << "Nodes: " << nodes / nFolds << std::endl; + // std::cout << "**********************************************************************************" << std::endl; + // std::cout << "Average Score Train: " << total_score_train / nFolds << std::endl; + // std::cout << "Average Score Test : " << total_score / nFolds << std::endl;return 0; } \ No newline at end of file diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index 5c8a536..717bf83 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -55,6 +55,7 @@ add_executable(b_main commands/b_main.cpp ${main_sources} common/Datasets.cpp common/Dataset.cpp common/Discretization.cpp reports/ReportConsole.cpp reports/ReportBase.cpp results/Result.cpp + experimental_clfs/XA1DE.cpp ) target_link_libraries(b_main "${PyClassifiers}" "${BayesNet}" fimdlp ${Python3_LIBRARIES} "${TORCH_LIBRARIES}" ${LIBTORCH_PYTHON} Boost::python Boost::numpy) diff --git a/src/commands/b_grid.cpp b/src/commands/b_grid.cpp index 2306b26..d93b689 100644 --- a/src/commands/b_grid.cpp +++ b/src/commands/b_grid.cpp @@ -8,7 +8,7 @@ #include "main/modelRegister.h" #include "main/ArgumentsExperiment.h" #include "common/Paths.h" -#include "common/Timer.h" +#include "common/Timer.hpp" #include "common/Colors.h" #include "common/DotEnv.h" #include "grid/GridSearch.h" diff --git a/src/common/Timer.h b/src/common/Timer.hpp similarity index 100% rename from src/common/Timer.h rename to src/common/Timer.hpp diff --git a/src/experimental_clfs/CountingSemaphore.hpp b/src/experimental_clfs/CountingSemaphore.hpp new file mode 100644 index 0000000..e217d00 --- /dev/null +++ b/src/experimental_clfs/CountingSemaphore.hpp @@ -0,0 +1,53 @@ +#ifndef COUNTING_SEMAPHORE_H +#define COUNTING_SEMAPHORE_H +#include +#include +#include +#include +#include +#include + +class CountingSemaphore { +public: + static CountingSemaphore& getInstance() + { + static CountingSemaphore instance; + return instance; + } + // Delete copy constructor and assignment operator + CountingSemaphore(const CountingSemaphore&) = delete; + CountingSemaphore& operator=(const CountingSemaphore&) = delete; + void acquire() + { + std::unique_lock lock(mtx_); + cv_.wait(lock, [this]() { return count_ > 0; }); + --count_; + } + void release() + { + std::lock_guard lock(mtx_); + ++count_; + if (count_ <= max_count_) { + cv_.notify_one(); + } + } + uint getCount() const + { + return count_; + } + uint getMaxCount() const + { + return max_count_; + } +private: + CountingSemaphore() + : max_count_(std::max(1u, static_cast(0.95 * std::thread::hardware_concurrency()))), + count_(max_count_) + { + } + std::mutex mtx_; + std::condition_variable cv_; + const uint max_count_; + uint count_; +}; +#endif \ No newline at end of file diff --git a/src/experimental_clfs/XA1DE.cpp b/src/experimental_clfs/XA1DE.cpp new file mode 100644 index 0000000..72d2ac2 --- /dev/null +++ b/src/experimental_clfs/XA1DE.cpp @@ -0,0 +1,150 @@ +// *************************************************************** +// SPDX-FileCopyrightText: Copyright 2025 Ricardo Montañana Gómez +// SPDX-FileType: SOURCE +// SPDX-License-Identifier: MIT +// *************************************************************** + +#include "XA1DE.h" + +namespace platform { + XA1DE::XA1DE() : semaphore_{ CountingSemaphore::getInstance() } + { + validHyperparameters = { "use_threads" }; + } + void XA1DE::setHyperparameters(const nlohmann::json& hyperparameters_) + { + auto hyperparameters = hyperparameters_; + if (hyperparameters.contains("use_threads")) { + use_threads = hyperparameters["use_threads"].get(); + hyperparameters.erase("use_threads"); + } + if (!hyperparameters.empty()) { + throw std::invalid_argument("Invalid hyperparameters" + hyperparameters.dump()); + } + } + void XA1DE::fit(std::vector> X, std::vector y, std::vector weights) + { + Timer timer, timert; + timer.start(); + timert.start(); + weights_ = weights; + std::vector> instances = X; + instances.push_back(y); + int num_instances = instances[0].size(); + int num_attributes = instances.size(); + normalize_weights(num_instances); + std::vector states; + for (int i = 0; i < num_attributes; i++) { + states.push_back(*max_element(instances[i].begin(), instances[i].end()) + 1); + } + aode_.init(states); + aode_.duration_first += timer.getDuration(); timer.start(); + std::vector instance; + for (int n_instance = 0; n_instance < num_instances; n_instance++) { + instance.clear(); + for (int feature = 0; feature < num_attributes; feature++) { + instance.push_back(instances[feature][n_instance]); + } + aode_.addSample(instance, weights_[n_instance]); + } + aode_.duration_second += timer.getDuration(); timer.start(); + // if (debug) aode_.show(); + aode_.computeProbabilities(); + aode_.duration_third += timer.getDuration(); + if (debug) { + // std::cout << "* Checking coherence... "; + // aode_.checkCoherenceApprox(1e-6); + // std::cout << "Ok!" << std::endl; + // aode_.show(); + // std::cout << "* Accumulated first time: " << aode_.duration_first << std::endl; + // std::cout << "* Accumulated second time: " << aode_.duration_second << std::endl; + // std::cout << "* Accumulated third time: " << aode_.duration_third << std::endl; + std::cout << "* Time to build the model: " << timert.getDuration() << " seconds" << std::endl; + // exit(1); + } + } + std::vector> XA1DE::predict_proba(std::vector>& test_data) + { + int test_size = test_data[0].size(); + std::vector> probabilities; + + std::vector instance; + for (int i = 0; i < test_size; i++) { + instance.clear(); + for (int j = 0; j < (int)test_data.size(); j++) { + instance.push_back(test_data[j][i]); + } + probabilities.push_back(aode_.predict_proba(instance)); + } + return probabilities; + } + std::vector> XA1DE::predict_proba_threads(const std::vector>& test_data) + { + int test_size = test_data[0].size(); + int sample_size = test_data.size(); + auto probabilities = std::vector>(test_size, std::vector(aode_.statesClass())); + + int chunk_size = std::min(150, int(test_size / semaphore_.getMaxCount()) + 1); + std::vector threads; + auto worker = [&](const std::vector>& samples, int begin, int chunk, int sample_size, std::vector>& predictions) { + std::string threadName = "(V)PWorker-" + std::to_string(begin) + "-" + std::to_string(chunk); +#if defined(__linux__) + pthread_setname_np(pthread_self(), threadName.c_str()); +#else + pthread_setname_np(threadName.c_str()); +#endif + + std::vector instance(sample_size); + for (int sample = begin; sample < begin + chunk; ++sample) { + for (int feature = 0; feature < sample_size; ++feature) { + instance[feature] = samples[feature][sample]; + } + predictions[sample] = aode_.predict_proba(instance); + } + semaphore_.release(); + }; + for (int begin = 0; begin < test_size; begin += chunk_size) { + int chunk = std::min(chunk_size, test_size - begin); + semaphore_.acquire(); + threads.emplace_back(worker, test_data, begin, chunk, sample_size, std::ref(probabilities)); + } + for (auto& thread : threads) { + thread.join(); + } + return probabilities; + } + std::vector XA1DE::predict(std::vector>& test_data) + { + auto probabilities = predict_proba(test_data); + std::vector predictions(probabilities.size(), 0); + + for (size_t i = 0; i < probabilities.size(); i++) { + predictions[i] = std::distance(probabilities[i].begin(), std::max_element(probabilities[i].begin(), probabilities[i].end())); + } + + return predictions; + } + float XA1DE::score(std::vector>& test_data, std::vector& labels) + { + aode_.duration_first = 0.0; + aode_.duration_second = 0.0; + aode_.duration_third = 0.0; + Timer timer; + timer.start(); + std::vector predictions = predict(test_data); + int correct = 0; + + for (size_t i = 0; i < predictions.size(); i++) { + if (predictions[i] == labels[i]) { + correct++; + } + } + if (debug) { + std::cout << "* Time to predict: " << timer.getDurationString() << std::endl; + std::cout << "* Accumulated first time: " << aode_.duration_first << std::endl; + std::cout << "* Accumulated second time: " << aode_.duration_second << std::endl; + std::cout << "* Accumulated third time: " << aode_.duration_third << std::endl; + } + return static_cast(correct) / predictions.size(); + } +} \ No newline at end of file diff --git a/src/experimental_clfs/XA1DE.h b/src/experimental_clfs/XA1DE.h new file mode 100644 index 0000000..6839f96 --- /dev/null +++ b/src/experimental_clfs/XA1DE.h @@ -0,0 +1,74 @@ +// *************************************************************** +// SPDX-FileCopyrightText: Copyright 2025 Ricardo Montañana Gómez +// SPDX-FileType: SOURCE +// SPDX-License-Identifier: MIT +// *************************************************************** + +#ifndef XA1DE_H +#define XA1DE_H +#include +#include +#include +#include +#include +#include "bayesnet/BaseClassifier.h" +#include "common/Timer.hpp" +#include "CountingSemaphore.hpp" +#include "Xaode.hpp" + +namespace platform { + class XA1DE : public bayesnet::BaseClassifier { + public: + XA1DE(); + virtual ~XA1DE() = default; + void setDebug(bool debug) { this->debug = debug; } + std::vector> predict_proba_threads(const std::vector>& test_data); + + XA1DE& fit(std::vector>& X, std::vector& y, const std::vector& features, const std::string& className, std::map>& states, const bayesnet::Smoothing_t smoothing) override; + XA1DE& fit(torch::Tensor& X, torch::Tensor& y, const std::vector& features, const std::string& className, std::map>& states, const bayesnet::Smoothing_t smoothing) override; + XA1DE& fit(torch::Tensor& dataset, const std::vector& features, const std::string& className, std::map>& states, const bayesnet::Smoothing_t smoothing) override; + XA1DE& fit(torch::Tensor& dataset, const std::vector& features, const std::string& className, std::map>& states, const torch::Tensor& weights, const bayesnet::Smoothing_t smoothing) override; + int getNumberOfNodes() const override { return 0; }; + int getNumberOfEdges() const override { return 0; }; + int getNumberOfStates() const override { return 0; }; + int getClassNumStates() const override { return 0; }; + torch::Tensor predict(torch::Tensor& X) override { return torch::zeros(0); }; + std::vector predict(std::vector>& X) override; + torch::Tensor predict_proba(torch::Tensor& X) override { return torch::zeros(0); }; + std::vector> predict_proba(std::vector>& X) override; + bayesnet::status_t getStatus() const override { return status; } + std::string getVersion() override { return { project_version.begin(), project_version.end() }; }; + float score(torch::Tensor& X, torch::Tensor& y) override { return 0; }; + float score(std::vector>& X, std::vector& y) override; + std::vector show() const override { return {}; } + std::vector topological_order() override { return {}; } + std::vector getNotes() const override { return notes; } + std::string dump_cpt() const override { return ""; } + void setHyperparameters(const nlohmann::json& hyperparameters) override; + + std::vector& getValidHyperparameters() { return validHyperparameters; } + protected: + void trainModel(const torch::Tensor& weights, const bayesnet::Smoothing_t smoothing) override; + + private: + inline void normalize_weights(int num_instances) + { + double sum = std::accumulate(weights_.begin(), weights_.end(), 0.0); + if (sum == 0) { + throw std::runtime_error("Weights sum zero."); + } + for (double& w : weights_) { + w = w * num_instances / sum; + } + } + // The instances of the dataset + Xaode aode_; + std::vector weights_; + CountingSemaphore& semaphore_; + bool debug = false; + bayesnet::status_t status = bayesnet::NORMAL; + std::vector notes; + bool use_threads = false; + }; +} +#endif // XA1DE_H \ No newline at end of file diff --git a/src/experimental_clfs/Xaode.hpp b/src/experimental_clfs/Xaode.hpp new file mode 100644 index 0000000..5093853 --- /dev/null +++ b/src/experimental_clfs/Xaode.hpp @@ -0,0 +1,579 @@ +// *************************************************************** +// SPDX-FileCopyrightText: Copyright 2025 Ricardo Montañana Gómez +// SPDX-FileType: SOURCE +// SPDX-License-Identifier: MIT +// *************************************************************** +// Based on the Geoff. I. Webb A1DE java algorithm +// https://weka.sourceforge.io/packageMetaData/AnDE/Latest.html + +#ifndef XAODE_H +#define XAODE_H +#include +#include +#include +#include +#include +#include +#include +#include + +namespace platform { + class Xaode { + public: + // ------------------------------------------------------- + // The Xaode can be EMPTY (just created), in COUNTS mode (accumulating raw counts) + // or PROBS mode (storing conditional probabilities). + enum class MatrixState { + EMPTY, + COUNTS, + PROBS + }; + double duration_first = 0.0; + double duration_second = 0.0; + double duration_third = 0.0; + Xaode() : nFeatures_{ 0 }, statesClass_{ 0 }, totalSize_{ 0 }, matrixState_{ MatrixState::EMPTY } {} + // ------------------------------------------------------- + // init + // ------------------------------------------------------- + // + // states.size() = nFeatures + 1, + // where states.back() = number of class states. + // + // We'll store: + // 1) p(c) in classPriors_ + // 2) p(x_i=si | c) in classFeatureProbs_ + // 3) p(x_j=sj | c, x_i=si) in data_, with i i is "superparent," j is "child." + // + // Internally, in COUNTS mode, data_ accumulates raw counts, then + // computeProbabilities(...) normalizes them into conditionals. + // + void init(const std::vector& states) + { + if (matrixState_ != MatrixState::EMPTY) { + throw std::logic_error("Xaode: already initialized."); + } + states_ = states; + nFeatures_ = static_cast(states_.size()) - 1; + if (nFeatures_ < 1) { + throw std::invalid_argument("Xaode: need at least 1 feature plus class states."); + } + statesClass_ = states_.back(); + if (statesClass_ <= 0) { + throw std::invalid_argument("Xaode: class states must be > 0."); + } + int totalStates = std::accumulate(states.begin(), states.end(), 0) - statesClass_; + + // For p(x_i=si | c), we store them in a 1D array classFeatureProbs_ after we compute. + // We'll need the offsets for each feature i in featureClassOffset_. + featureClassOffset_.resize(nFeatures_); + // We'll store p(x_child=sj | c, x_sp=si) for each pair (i& instance, double weight) + { + // + // (A) increment classCounts_ + // (B) increment feature–class counts => for p(x_i|c) + // (C) increment pair (superparent= i, child= j) counts => data_ + // + + // if (matrixState_ != MatrixState::COUNTS) { + // throw std::logic_error("addSample: not in COUNTS mode."); + // } + // if (static_cast(instance.size()) != nFeatures_ + 1) { + // throw std::invalid_argument("addSample: instance.size() must be nFeatures_ + 1."); + // } + + int c = instance.back(); + // if (c < 0 || c >= statesClass_) { + // throw std::out_of_range("addSample: class index out of range."); + // } + if (weight <= 0.0) { + return; + } + // (A) increment classCounts_ + classCounts_[c] += weight; + + // (B,C) + // We'll store raw counts now and turn them into p(child| c, superparent) later. + int idx, fcIndex, si, sj, i_offset; + for (int i = 0; i < nFeatures_; ++i) { + si = instance[i]; + // (B) increment feature–class counts => for p(x_i|c) + fcIndex = (featureClassOffset_[i] + si) * statesClass_ + c; + classFeatureCounts_[fcIndex] += weight; + // (C) increment pair (superparent= i, child= j) counts => data_ + i_offset = pairOffset_[featureClassOffset_[i] + si]; + for (int j = 0; j < i; ++j) { + sj = instance[j]; + idx = (i_offset + featureClassOffset_[j] + sj) * statesClass_ + c; + data_[idx] += weight; + } + } + } + + // ------------------------------------------------------- + // computeProbabilities + // ------------------------------------------------------- + // + // Once all samples are added in COUNTS mode, call this to: + // 1) compute class priors p(c) + // 2) compute p(x_i=si | c) => classFeatureProbs_ + // 3) compute p(x_j=sj | c, x_i=si) => data_ (for ij) + // + void computeProbabilities() + { + if (matrixState_ != MatrixState::COUNTS) { + throw std::logic_error("computeProbabilities: must be in COUNTS mode."); + } + // (1) p(c) + double totalCount = std::accumulate(classCounts_.begin(), classCounts_.end(), 0.0); + if (totalCount <= 0.0) { + // fallback => uniform + double unif = 1.0 / statesClass_; + for (int c = 0; c < statesClass_; ++c) { + classPriors_[c] = unif; + } + } else { + for (int c = 0; c < statesClass_; ++c) { + classPriors_[c] = classCounts_[c] / totalCount; + } + } + // (2) p(x_i=si | c) => classFeatureProbs_ + int idx, sf; + double denom, countVal, p; + for (int feature = 0; feature < nFeatures_; ++feature) { + sf = states_[feature]; + for (int c = 0; c < statesClass_; ++c) { + denom = classCounts_[c] * sf; + if (denom <= 0.0) { + // fallback => uniform + for (int sf_value = 0; sf_value < sf; ++sf_value) { + idx = (featureClassOffset_[feature] + sf_value) * statesClass_ + c; + classFeatureProbs_[idx] = 1.0 / sf; + } + } else { + for (int sf_value = 0; sf_value < sf; ++sf_value) { + idx = (featureClassOffset_[feature] + sf_value) * statesClass_ + c; + countVal = classFeatureCounts_[idx]; + p = ((countVal + SMOOTHING / (statesClass_ * states_[feature])) / (totalCount + SMOOTHING)); + classFeatureProbs_[idx] = p; + } + } + } + } + // getCountFromTable(int classVal, int pIndex, int childIndex) + // (3) p(x_j=sj | c, x_i=si) => data_(i,si,j,sj,c) + // (3) p(x_i=si | c, x_j=sj) => dataOpp_(j,sj,i,si,c) + double pccCount, pcCount, ccCount; + double conditionalProb, oppositeCondProb; + int part1, part2, p1, part2_class, p1_class; + for (int parent = nFeatures_ - 1; parent >= 0; --parent) { + // for (int parent = 3; parent >= 3; --parent) { + for (int sp = 0; sp < states_[parent]; ++sp) { + p1 = featureClassOffset_[parent] + sp; + part1 = pairOffset_[p1]; + p1_class = p1 * statesClass_; + for (int child = parent - 1; child >= 0; --child) { + // for (int child = 2; child >= 2; --child) { + for (int sc = 0; sc < states_[child]; ++sc) { + part2 = featureClassOffset_[child] + sc; + part2_class = part2 * statesClass_; + for (int c = 0; c < statesClass_; c++) { + //idx = compute_index(parent, sp, child, sc, classval); + idx = (part1 + part2) * statesClass_ + c; + // Parent, Child, Class Count + pccCount = data_[idx]; + // Parent, Class count + pcCount = classFeatureCounts_[p1_class + c]; + // Child, Class count + ccCount = classFeatureCounts_[part2_class + c]; + conditionalProb = (pccCount + SMOOTHING / states_[parent]) / (ccCount + SMOOTHING); + data_[idx] = conditionalProb; + oppositeCondProb = (pccCount + SMOOTHING / states_[child]) / (pcCount + SMOOTHING); + dataOpp_[idx] = oppositeCondProb; + } + } + } + } + } + matrixState_ = MatrixState::PROBS; + } + + // ------------------------------------------------------- + // predict_proba_spode + // ------------------------------------------------------- + // + // Single-superparent approach: + // P(c | x) ∝ p(c) * p(x_sp| c) * ∏_{i≠sp} p(x_i | c, x_sp) + // + // 'instance' should have size == nFeatures_ (no class). + // sp in [0..nFeatures_). + // We multiply p(c) * p(x_sp| c) * p(x_i| c, x_sp). + // Then normalize the distribution. + // + std::vector predict_proba_spode(const std::vector& instance, int parent) const + { + if (matrixState_ != MatrixState::PROBS) { + throw std::logic_error("predict_proba_spode: Xaode not in PROBS state."); + } + if ((int)instance.size() != nFeatures_) { + throw std::invalid_argument("predict_proba_spode: instance.size() != nFeatures_."); + } + if (parent < 0 || parent >= nFeatures_) { + throw std::out_of_range("predict_proba_spode: invalid superparent index."); + } + + std::vector scores(statesClass_, 0.0); + int sp = instance[parent]; + int idx; + double pSpGivenC, pChildGivenSp, product; + double base; + double offset = (featureClassOffset_[parent] + sp) * statesClass_; + double parent_offset = pairOffset_[featureClassOffset_[parent] + sp]; + // For each class c + for (int c = 0; c < statesClass_; ++c) { + // Start with p(c) * p(x_sp=spState| c) + pSpGivenC = classFeatureProbs_[offset + c]; + product = pSpGivenC; + bool zeroProb = false; + for (int feature = 0; feature < nFeatures_; ++feature) { + if (feature == parent) continue; + int sf = instance[feature]; + // Retrieve p(x_i= state_i | c, x_sp= spState) + base = (parent_offset + featureClassOffset_[feature] + sf) * statesClass_; + idx = base + c; + pChildGivenSp = data_[idx] * dataOpp_[idx]; + if (pChildGivenSp <= 0.0) { + zeroProb = true; + break; + } + product *= pChildGivenSp; + } + scores[c] = zeroProb ? 0.0 : product; + } + normalize(scores); + return scores; + } + std::vector predict_proba(std::vector& instance) + { + Timer timer; + timer.start(); + if (matrixState_ != MatrixState::PROBS) { + throw std::logic_error("predict_proba: Xaode not in PROBS state."); + } + if ((int)instance.size() != nFeatures_) { + throw std::invalid_argument("predict_proba: instance.size() != nFeatures_."); + } + // accumulates posterior probabilities for each class + auto probs = std::vector(statesClass_); + auto spodeProbs = std::vector>(nFeatures_, std::vector(statesClass_)); + // Initialize the probabilities with the feature|class probabilities + int localOffset; + for (int feature = 0; feature < nFeatures_; ++feature) { + localOffset = (featureClassOffset_[feature] + instance[feature]) * statesClass_; + for (int c = 0; c < statesClass_; ++c) { + spodeProbs[feature][c] = classFeatureProbs_[localOffset + c]; + } + } + duration_first += timer.getDuration(); timer.start(); + int idx, base, sp, sc, parent_offset; + for (int parent = 1; parent < nFeatures_; ++parent) { + sp = instance[parent]; + parent_offset = pairOffset_[featureClassOffset_[parent] + sp]; + for (int child = 0; child < parent; ++child) { + sc = instance[child]; + base = (parent_offset + featureClassOffset_[child] + sc) * statesClass_; + for (int c = 0; c < statesClass_; ++c) { + /* + * The probability P(xc|xp,c) is stored in dataOpp_, and + * the probability P(xp|xc,c) is stored in data_ + */ + /* + int base = pairOffset_[i * nFeatures_ + j]; + int blockSize = states_[i] * states_[j]; + return base + c * blockSize + (si * states_[j] + sj); + */ + // index = compute_index(parent, instance[parent], child, instance[child], classVal); + idx = base + c; + spodeProbs[child][c] *= data_[idx]; + // spodeProbs[child][c] *= data_.at(index); + spodeProbs[parent][c] *= dataOpp_[idx]; + // spodeProbs[parent][c] *= dataOpp_.at(index); + } + } + } + duration_second += timer.getDuration(); timer.start(); + /* add all the probabilities for each class */ + for (int c = 0; c < statesClass_; ++c) { + for (int i = 0; i < nFeatures_; ++i) { + probs[c] += spodeProbs[i][c]; + } + } + // Normalize the probabilities + normalize(probs); + return probs; + } + void normalize(std::vector& probs) const + { + double sum = 0; + for (double d : probs) { + sum += d; + } + if (std::isnan(sum)) { + throw std::runtime_error("Can't normalize array. Sum is NaN."); + } + if (sum == 0) { + return; + } + for (int i = 0; i < (int)probs.size(); i++) { + probs[i] /= sum; + } + } + + // ------------------------------------------------------- + // checkCoherence + // ------------------------------------------------------- + // + // Check that the class priors, feature–class distributions and pairwise conditionals + // are coherent. They have to sum to 1.0 within a threshold. + // + void checkCoherenceApprox(double threshold) const + { + if (matrixState_ != MatrixState::PROBS) { + throw std::logic_error("checkCoherenceApprox: must be in PROBS state."); + } + + // ------------------------------------------------------------------ + // 1) Check that sum of class priors ~ 1 + // ------------------------------------------------------------------ + double sumPriors = 0.0; + for (double pc : classPriors_) { + sumPriors += pc; + } + if (std::fabs(sumPriors - 1.0) > threshold) { + std::ostringstream oss; + oss << "Xaode::checkCoherenceApprox - sum of classPriors = " << sumPriors + << ", differs from 1.0 by more than " << threshold; + throw std::runtime_error(oss.str()); + } + + // ------------------------------------------------------------------ + // 2) For each feature i and class c, the sum over all states si of + // classFeatureProbs_ should match the prior p(c) ~ classPriors_[c]. + // + // (Because if you're storing p(x_i=si, c)/total or a scaled version, + // summing over si is effectively p(c).) + // ------------------------------------------------------------------ + for (int c = 0; c < statesClass_; ++c) { + for (int i = 0; i < nFeatures_; ++i) { + double sumFeature = 0.0; + for (int si = 0; si < states_[i]; ++si) { + int idx = (featureClassOffset_[i] + si) * statesClass_ + c; + sumFeature += classFeatureProbs_[idx]; + } + double expected = classPriors_[c]; + if (std::fabs(sumFeature - expected) > threshold) { + std::ostringstream oss; + oss << "Xaode::checkCoherenceApprox - sum_{si} classFeatureProbs_ " + << "for (feature=" << i << ", class=" << c << ") = " << sumFeature + << ", expected ~ " << expected + << ", difference is " << std::fabs(sumFeature - expected) + << " > threshold=" << threshold; + throw std::runtime_error(oss.str()); + } + } + } + + // ------------------------------------------------------------------ + // 3) For data_: sum_{child states} data_ should match the "parent" row + // in classFeatureProbs_, i.e. p(x_i=si, c). + // + // Because if data_[... i, si, j, sj, c] holds something like + // p(x_i=si, x_j=sj, c) (or a scaled fraction), + // then sum_{ sj } data_ = p(x_i=si, c). + // ------------------------------------------------------------------ + for (int parent = 1; parent < nFeatures_; ++parent) { + for (int child = 0; child < parent; ++child) { + for (int c = 0; c < statesClass_; ++c) { + for (int spVal = 0; spVal < states_[parent]; ++spVal) { + double sumChildProb = 0.0; + // pairOffset_ gives the offset for (parent featureVal), + // then we add the child's offset and multiply by statesClass_. + int part1 = pairOffset_[featureClassOffset_[parent] + spVal]; + for (int scVal = 0; scVal < states_[child]; ++scVal) { + int part2 = featureClassOffset_[child] + scVal; + int idx = (part1 + part2) * statesClass_ + c; + sumChildProb += data_[idx]; + } + // Compare with classFeatureProbs_[parent, spVal, c] + double expected = classFeatureProbs_[ + (featureClassOffset_[parent] + spVal) * statesClass_ + c + ]; + if (std::fabs(sumChildProb - expected) > threshold) { + std::ostringstream oss; + oss << "Xaode::checkCoherenceApprox - sum_{sj} data_ " + << "for (parentFeature=" << parent + << ", parentVal=" << spVal + << ", childFeature=" << child + << ", class=" << c << ") = " << sumChildProb + << ", expected ~ " << expected + << ", diff " << std::fabs(sumChildProb - expected) + << " > threshold=" << threshold; + throw std::runtime_error(oss.str()); + } + } + } + } + } + + // ------------------------------------------------------------------ + // 4) For dataOpp_: sum_{parent states} dataOpp_ should match the "child" + // row in classFeatureProbs_, i.e. p(x_j=sj, c). + // ------------------------------------------------------------------ + for (int parent = 1; parent < nFeatures_; ++parent) { + for (int child = 0; child < parent; ++child) { + for (int c = 0; c < statesClass_; ++c) { + for (int scVal = 0; scVal < states_[child]; ++scVal) { + double sumParentProb = 0.0; + int part2 = featureClassOffset_[child] + scVal; + for (int spVal = 0; spVal < states_[parent]; ++spVal) { + int part1 = pairOffset_[featureClassOffset_[parent] + spVal]; + int idx = (part1 + part2) * statesClass_ + c; + sumParentProb += dataOpp_[idx]; + } + // Compare with classFeatureProbs_[child, scVal, c] + double expected = classFeatureProbs_[ + (featureClassOffset_[child] + scVal) * statesClass_ + c + ]; + if (std::fabs(sumParentProb - expected) > threshold) { + std::ostringstream oss; + oss << "Xaode::checkCoherenceApprox - sum_{spVal} dataOpp_ " + << "for (childFeature=" << child + << ", childVal=" << scVal + << ", parentFeature=" << parent + << ", class=" << c << ") = " << sumParentProb + << ", expected ~ " << expected + << ", diff " << std::fabs(sumParentProb - expected) + << " > threshold=" << threshold; + throw std::runtime_error(oss.str()); + } + } + } + } + } + + // If we get here, all sums are coherent under this "joint distribution" interpretation + } + int statesClass() const + { + return statesClass_; + } + + private: + // ----------- + // MEMBER DATA + // ----------- + std::vector states_; // [states_feat0, ..., states_feat(n-1), statesClass_] + int nFeatures_; + int statesClass_; + + // data_ means p(child=sj | c, superparent= si) after normalization. + // But in COUNTS mode, it accumulates raw counts. + std::vector pairOffset_; + int totalSize_; + // data_ stores p(child=sj | c, superparent=si) for each pair (i data_; + // dataOpp_ stores p(superparent=si | c, child=sj) for each pair (i dataOpp_; + + // classCounts_[c] + std::vector classCounts_; + + // For p(x_i=si| c), we store counts in classFeatureCounts_ => offset by featureClassOffset_[i] + std::vector featureClassOffset_; + std::vector classFeatureCounts_; + std::vector classFeatureProbs_; // => p(x_i=si | c) after normalization + + std::vector classPriors_; // => p(c) + + MatrixState matrixState_; + + double SMOOTHING = 1.0; + }; +} +#endif // XAODE_H \ No newline at end of file diff --git a/src/grid/GridBase.h b/src/grid/GridBase.h index e496bca..8ae3421 100644 --- a/src/grid/GridBase.h +++ b/src/grid/GridBase.h @@ -5,7 +5,7 @@ #include #include #include "common/Datasets.h" -#include "common/Timer.h" +#include "common/Timer.hpp" #include "common/Colors.h" #include "main/HyperParameters.h" #include "GridData.h" diff --git a/src/grid/GridConfig.h b/src/grid/GridConfig.h index a9159f0..6a1acb5 100644 --- a/src/grid/GridConfig.h +++ b/src/grid/GridConfig.h @@ -5,7 +5,7 @@ #include #include #include "common/Datasets.h" -#include "common/Timer.h" +#include "common/Timer.hpp" #include "main/HyperParameters.h" #include "GridData.h" #include "GridConfig.h" diff --git a/src/grid/GridSearch.h b/src/grid/GridSearch.h index 6f8ab37..2797b78 100644 --- a/src/grid/GridSearch.h +++ b/src/grid/GridSearch.h @@ -6,7 +6,7 @@ #include #include #include "common/Datasets.h" -#include "common/Timer.h" +#include "common/Timer.hpp" #include "main/HyperParameters.h" #include "GridData.h" #include "GridBase.h" diff --git a/src/main/Models.h b/src/main/Models.h index e8eb58f..7e8e1b4 100644 --- a/src/main/Models.h +++ b/src/main/Models.h @@ -20,6 +20,7 @@ #include #include #include +#include "../experimental_clfs/XA1DE.h" namespace platform { class Models { public: diff --git a/src/reports/ReportConsole.cpp b/src/reports/ReportConsole.cpp index ef0eb97..55037ce 100644 --- a/src/reports/ReportConsole.cpp +++ b/src/reports/ReportConsole.cpp @@ -2,7 +2,7 @@ #include #include "best/BestScore.h" #include "common/CLocale.h" -#include "common/Timer.h" +#include "common/Timer.hpp" #include "ReportConsole.h" #include "main/Scores.h" @@ -251,7 +251,7 @@ namespace platform { if (train_data) { oss << color_line << std::left << std::setw(maxLine) << output_train[i] << suffix << Colors::BLUE() << " | " << color_line << std::left << std::setw(maxLine) - << output_test[i] << std::endl; + << output_test[i] << std::endl; } else { oss << color_line << output_test[i] << std::endl; } diff --git a/src/results/Result.h b/src/results/Result.h index 78bfd62..3f45c70 100644 --- a/src/results/Result.h +++ b/src/results/Result.h @@ -4,7 +4,7 @@ #include #include #include -#include "common/Timer.h" +#include "common/Timer.hpp" #include "main/HyperParameters.h" #include "main/PartialResult.h"