Skip to content

Commit

Permalink
started creating more classes
Browse files Browse the repository at this point in the history
  • Loading branch information
osorensen committed Nov 21, 2023
1 parent 7840bd0 commit 7c650d7
Show file tree
Hide file tree
Showing 5 changed files with 137 additions and 86 deletions.
2 changes: 1 addition & 1 deletion src/mixtures.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,7 @@ uvec update_cluster_labels(
const mat& dist_mat,
const vec& cluster_probs,
const vec& alpha_old,
const int& n_items,
const unsigned int& n_items,
const int& t,
const std::string& metric,
const Rcpp::List& logz_list,
Expand Down
2 changes: 1 addition & 1 deletion src/mixtures.h
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,7 @@ arma::uvec update_cluster_labels(
const arma::mat& dist_mat,
const arma::vec& cluster_probs,
const arma::vec& alpha_old,
const int& n_items,
const unsigned int& n_items,
const int& t,
const std::string& metric,
const Rcpp::List& logz_list,
Expand Down
63 changes: 39 additions & 24 deletions src/parameters.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -2,43 +2,54 @@

using namespace arma;

Data::Data(
const Rcpp::List& data
) :
rankings { Rcpp::as<mat>(data["rankings"]).t() },
n_assessors { rankings.n_cols },
n_items { rankings.n_rows }
{

}

Priors::Priors(
const Rcpp::List& priors
) : lambda { verify_positive(Rcpp::as<double>(priors["lambda"])) },
kappa_1 { Rcpp::as<unsigned int>(priors["kappa_1"]) },
kappa_2 { Rcpp::as<unsigned int>(priors["kappa_2"]) }
{

}

Parameters::Parameters(
const Rcpp::List& model,
const Rcpp::List& compute_options,
const Rcpp::List& priors,
const Rcpp::List& initial_values,
const int n_items) :
metric { verify_metric(model["metric"]) },
n_items { n_items },
nmc { verify_positive(compute_options["nmc"]) }
const unsigned int n_items) :
alpha_jump { Rcpp::as<int>(compute_options["alpha_jump"]) },
alpha_prop_sd { verify_positive(Rcpp::as<double>(compute_options["alpha_prop_sd"])) },
error_model { verify_error_model(Rcpp::as<std::string>(model["error_model"])) },
leap_size { Rcpp::as<int>(compute_options["leap_size"]) },
metric { verify_metric(Rcpp::as<std::string>(model["metric"])) },
n_clusters { Rcpp::as<int>(model["n_clusters"]) },
nmc { Rcpp::as<int>(compute_options["nmc"]) },
rho_thinning { Rcpp::as<int>(compute_options["rho_thinning"]) }
{
int n_clusters = model["n_clusters"];
int alpha_jump = compute_options["alpha_jump"];

alpha.set_size(n_clusters, std::ceil(static_cast<double>(nmc * 1.0 / alpha_jump)));
double alpha_init = initial_values["alpha_init"];
alpha.col(0).fill(alpha_init);
alpha_old = alpha.col(0);
lambda = Rcpp::as<double>(priors["lambda"]);
alpha_prop_sd = Rcpp::as<double>(compute_options["alpha_prop_sd"]);

rho_thinning = Rcpp::as<int>(compute_options["rho_thinning"]);
rho.set_size(n_items, n_clusters, std::ceil(static_cast<double>(nmc * 1.0 / rho_thinning)));
Rcpp::Nullable<mat> rho_init = initial_values["rho_init"];
rho.slice(0) = initialize_rho(n_items, n_clusters, rho_init);
rho_old = rho(span::all, span::all, span(0));

leap_size = Rcpp::as<int>(compute_options["leap_size"]);

kappa_1 = Rcpp::as<int>(priors["kappa_1"]);
kappa_2 = Rcpp::as<int>(priors["kappa_2"]);

error_model = Rcpp::as<std::string>(model["error_model"]);
if(error_model == "bernoulli"){
theta = zeros<vec>(nmc);
shape_1 = zeros<vec>(nmc);
shape_2 = zeros<vec>(nmc);
shape_1(0) = kappa_1;
shape_2(0) = kappa_2;
} else {
theta.reset();
shape_1.reset();
Expand All @@ -61,9 +72,10 @@ void Parameters::update_rho(int cluster_index, int t, int& rho_index,
}

void Parameters::update_shape(int t, const mat& rankings,
const Rcpp::List& constraints) {

const Rcpp::List& constraints,
const Priors& priors) {

const unsigned int n_items = rankings.n_rows;
int n_assessors = rankings.n_cols;
int sum_1{};
int sum_2{};
Expand All @@ -85,8 +97,9 @@ void Parameters::update_shape(int t, const mat& rankings,
}
}
}
shape_1(t) = kappa_1 + sum_1;
shape_2(t) = kappa_2 + sum_2;

shape_1(t) = priors.kappa_1 + sum_1;
shape_2(t) = priors.kappa_2 + sum_2;
theta(t) = rtruncbeta(shape_1(t), shape_2(t), 0.5);
}

Expand All @@ -95,8 +108,10 @@ void Parameters::update_alpha(
int alpha_index,
const mat& rankings,
const vec& observation_frequency,
const Rcpp::List& logz_list) {
const Rcpp::List& logz_list,
const Priors& priors) {

const unsigned int n_items = rankings.n_rows;
double alpha_proposal = std::exp(randn<double>() * alpha_prop_sd +
std::log(alpha_old(cluster_index)));

Expand All @@ -108,7 +123,7 @@ void Parameters::update_alpha(
// Compute the Metropolis-Hastings ratio
double ratio =
alpha_diff / n_items * rank_dist +
lambda * alpha_diff +
priors.lambda * alpha_diff +
sum(observation_frequency) * (
get_partition_function(n_items, alpha_old(cluster_index), logz_list, metric) -
get_partition_function(n_items, alpha_proposal, logz_list, metric)
Expand Down
83 changes: 62 additions & 21 deletions src/parameters.h
Original file line number Diff line number Diff line change
Expand Up @@ -7,16 +7,43 @@
#include "distances.h"
#include "partitionfuns.h"

template <typename T>
static T verify_positive(const T input) {
if(input < 0) {
Rcpp::stop("Positive value required.\n");
}
return input;
}

struct Data {
Data(const Rcpp::List& data);
~Data() = default;

arma::mat rankings;
const unsigned int n_assessors;
const unsigned int n_items;

};

struct Priors {
Priors(const Rcpp::List& priors);
~Priors() = default;

const double lambda;
const unsigned int kappa_1;
const unsigned int kappa_2;
};

struct Parameters {
Parameters(const Rcpp::List& model,
const Rcpp::List& compute_options,
const Rcpp::List& priors,
const Rcpp::List& initial_values,
const int n_items);
const unsigned int n_items);
~Parameters() = default;

void update_shape(int t, const arma::mat& rankings,
const Rcpp::List& constraints);
const Rcpp::List& constraints,
const Priors& priors);
void update_rho(int cluster_index, int t, int& rho_index,
const arma::mat& rankings,
const arma::vec& observation_frequency);
Expand All @@ -26,31 +53,35 @@ struct Parameters {
int alpha_index,
const arma::mat& rankings,
const arma::vec& observation_frequency,
const Rcpp::List& logz_list);
const Rcpp::List& logz_list,
const Priors& priors);

arma::mat alpha;
arma::vec alpha_old;
arma::cube rho;
arma::mat rho_old;
std::string error_model;
arma::vec theta;
arma::vec shape_1;
arma::vec shape_2;
arma::vec theta;

int get_nmc() {
const int get_alpha_jump() {
return alpha_jump;
}
const int get_nmc() {
return nmc;
}
std::string get_metric() {
const int get_n_clusters() {
return n_clusters;
}
const std::string get_error_model() {
return error_model;
}
const std::string get_metric() {
return metric;
}

private:
static int verify_positive(const int input) {
if(input < 0) {
Rcpp::stop("Positive value required.\n");
}
return input;
}

static std::string verify_metric(const std::string input) {
bool check = (input.compare("footrule") == 0) ||
(input.compare("spearman") == 0) ||
Expand All @@ -63,15 +94,25 @@ struct Parameters {
}
return input;
}
int kappa_1;
int kappa_2;

static std::string verify_error_model(const std::string input) {
bool check = (input.compare("none") == 0) ||
(input.compare("bernoulli") == 0);
if(!check) {
Rcpp::stop("Unknown error model.\n");
}
return input;
}

const int alpha_jump;
const double alpha_prop_sd;
const std::string error_model;
const int leap_size;
const std::string metric;
const int n_items;
const int n_clusters;
const int nmc;
int leap_size;
int rho_thinning;
double lambda;
double alpha_prop_sd;
const int rho_thinning;

};


Expand Down
Loading

0 comments on commit 7c650d7

Please sign in to comment.