From 01f4f2fa3f72288ca89ffc142e7050ec2934ec2c Mon Sep 17 00:00:00 2001 From: Dongryeol Lee Date: Sun, 14 Oct 2007 18:11:27 +0000 Subject: [PATCH] Miscellaenous cleanup --- fastlib/u/dongryel/kde/kde.h | 119 ++++++++++++++++++++++++--------- fastlib/u/dongryel/kde/main.cc | 17 +++-- 2 files changed, 100 insertions(+), 36 deletions(-) diff --git a/fastlib/u/dongryel/kde/kde.h b/fastlib/u/dongryel/kde/kde.h index 69020ae07c..4d3d107d18 100644 --- a/fastlib/u/dongryel/kde/kde.h +++ b/fastlib/u/dongryel/kde/kde.h @@ -13,10 +13,14 @@ class NaiveKde { /** query dataset */ Matrix qset_; - + + ArrayList q_old_from_new_; + /** reference dataset */ Matrix rset_; + ArrayList r_old_from_new_; + /** kernel */ TKernel kernel_; @@ -56,18 +60,21 @@ class NaiveKde { densities_.SetZero(); } - void Init(Matrix &qset, Matrix &rset) { + void Init(Matrix &qset, ArrayList &q_old_from_new, + Matrix &rset, ArrayList &r_old_from_new) { // get datasets qset_.Alias(qset); + q_old_from_new_.Steal(&q_old_from_new); rset_.Alias(rset); + r_old_from_new_.Steal(&r_old_from_new); // get bandwidth kernel_.Init(fx_param_double_req(NULL, "bandwidth")); // allocate density storage densities_.Init(qset.n_cols()); - densities_.SetZero(); + densities_.SetZero(); } void PrintDebug() { @@ -78,7 +85,9 @@ class NaiveKde { if((fname = fx_param_str(NULL, "naive_kde_output", NULL)) != NULL) { stream = fopen(fname, "w+"); } - densities_.PrintDebug("Naive KDE results", stream); + for(index_t q = 0; q < qset_.n_cols(); q++) { + fprintf(stream, "%g\n", densities_[q_old_from_new_[q]]); + } if(stream != stdout) { fclose(stream); @@ -197,8 +206,11 @@ class FastKde { } void MergeChildBounds(KdeStat &left_stat, KdeStat &right_stat) { + + // steal left and right children's tokens double min_mass_t = min(left_stat.mass_t_, right_stat.mass_t_); + // improve lower and upper bound mass_l_ = max(mass_l_, min(left_stat.mass_l_, right_stat.mass_l_)); mass_u_ = min(mass_u_, max(left_stat.mass_u_, right_stat.mass_u_)); mass_t_ += min_mass_t; @@ -244,15 +256,27 @@ class FastKde { /** query dataset */ Matrix qset_; + /** + * maps new indices of the query dataset to the old indices (of original + * text dataset + */ + ArrayList q_old_from_new_; + /** query tree */ Tree *qroot_; /** reference dataset */ Matrix rset_; + /** + * maps new indices of the reference dataset to the old indices (of original + * text dataset + */ + ArrayList r_old_from_new_; + /** reference tree */ Tree *rroot_; - + /** reference weights */ Vector rset_weights_; @@ -288,7 +312,6 @@ class FastKde { // member functions void UpdateBounds(Tree *qnode, Tree *rnode, - KdeStat *left_stat, KdeStat *right_stat, double *dl, double *de, double *du, double *dt) { // query self statistics @@ -348,7 +371,7 @@ class FastKde { // for a leaf node, incorporate the lower and upper bound changes into // its additional offset - if(left_stat == NULL) { + if(qnode->is_leaf()) { qstat.more_l_ += dl_ref; qstat.more_u_ += du_ref; } @@ -357,10 +380,10 @@ class FastKde { // the immediate descendants else { - left_stat->owed_l_ += dl_ref; - left_stat->owed_u_ += du_ref; - right_stat->owed_l_ += dl_ref; - right_stat->owed_u_ += du_ref; + qnode->left()->stat().owed_l_ += dl_ref; + qnode->left()->stat().owed_u_ += du_ref; + qnode->right()->stat().owed_l_ += dl_ref; + qnode->right()->stat().owed_u_ += du_ref; } // zero out lower and upper bounds @@ -375,9 +398,14 @@ class FastKde { // compute unnormalized sum for(index_t q = qnode->begin(); q < qnode->end(); q++) { + // get query point const double *q_col = qset_.GetColumnPtr(q); for(index_t r = rnode->begin(); r < rnode->end(); r++) { + + // get reference point const double *r_col = rset_.GetColumnPtr(r); + + // pairwise distance and kernel value double dsqd = la::DistanceSqEuclidean(qset_.n_rows(), q_col, r_col); double ker_value = kernel_.EvalUnnormOnSq(dsqd); @@ -387,7 +415,7 @@ class FastKde { } } - // tally up the unused error components + // tally up the unused error components due to exhaustive computation qnode->stat().mass_t_ += rnode->count(); // get a tighter lower and upper bound by looping over each query point @@ -525,7 +553,7 @@ class FastKde { return 1; } else { - dl_ = du_ = dt_ = 0; + dl_ = de_ = du_ = dt_ = 0; return 0; } } @@ -550,31 +578,33 @@ class FastKde { /** canonical fast KDE case */ void FKde(Tree *qnode, Tree *rnode) { + // query node statistics KdeStat &stat = qnode->stat(); + + // left child and right child of query node statistics KdeStat *left_stat = NULL; KdeStat *right_stat = NULL; - // process density bound changes sent from the ancestor query nodes - // then tighten lower/upper bounds and the error reclaimed based on - // the children + // process density bound changes sent from the ancestor query nodes, + UpdateBounds(qnode, rnode, &stat.owed_l_, NULL, &stat.owed_u_, NULL); + + // for non-leaf query node, tighten lower/upper bounds and the + // reclaim tokens unused by the children. if(!qnode->is_leaf()) { left_stat = &(qnode->left()->stat()); right_stat = &(qnode->right()->stat()); - UpdateBounds(qnode, rnode, left_stat, right_stat, &stat.owed_l_, NULL, - &stat.owed_u_, NULL); stat.MergeChildBounds(*left_stat, *right_stat); } + // try finite difference pruning first if(Prunable(qnode, rnode)) { - UpdateBounds(qnode, rnode, left_stat, right_stat, &dl_, &de_, &du_, - &dt_); + UpdateBounds(qnode, rnode, &dl_, &de_, &du_, &dt_); return; } - // if prunable, then prune + // try series-expansion pruning else if(PrunableEnhanced(qnode, rnode)) { - UpdateBounds(qnode, rnode, left_stat, right_stat, &dl_, NULL, &du_, - &dt_); + UpdateBounds(qnode, rnode, &dl_, NULL, &du_, &dt_); return; } @@ -641,20 +671,35 @@ class FastKde { */ void PreProcess(Tree *node) { + // initialize the center of expansions and bandwidth for + // series expansion node->stat().Init(sqrt(kernel_.bandwidth_sq()), &sea_); node->bound().CalculateMidpoint (&(node->stat().farfield_expansion_.get_center())); node->bound().CalculateMidpoint (&(node->stat().local_expansion_.get_center())); + // initialize lower bound to 0 + node->stat().mass_l_ = 0; + + // set the finite difference approximated amounts to 0 + node->stat().mass_e_ = 0; + // set the upper bound to the number of reference points node->stat().mass_u_ = rset_.n_cols(); + // set the number of tokens to 0 + node->stat().mass_t_ = 0; + // for non-leaf node, recurse if(!node->is_leaf()) { + node->stat().owed_l_ = node->stat().owed_u_ = 0; PreProcess(node->left()); PreProcess(node->right()); } + else { + node->stat().more_l_ = node->stat().more_u_ = 0; + } } /** post processing step */ @@ -711,10 +756,18 @@ class FastKde { /** get the reference dataset */ Matrix &get_reference_dataset() { return rset_; } - + + ArrayList &get_reference_old_from_new_mapping() { + return r_old_from_new_; + } + /** get the query dataset */ Matrix &get_query_dataset() { return qset_; } + ArrayList &get_query_old_from_new_mapping() { + return q_old_from_new_; + } + /** get the density estimate */ const Vector &get_density_estimates() { return densities_e_; } @@ -728,7 +781,7 @@ class FastKde { // initialize the lower and upper bound densities densities_l_.SetZero(); densities_e_.SetZero(); - densities_u_.SetZero(); + densities_u_.SetAll(rset_.n_cols()); printf("\nStarting fast KDE...\n"); fx_timer_start(NULL, "fast_kde_compute"); @@ -744,7 +797,7 @@ class FastKde { order_farfield_ = order_local_ = order_farfield_to_local_ = -1; FKde(qroot_, rroot_); - // postprocessing step + // postprocessing step for finalizing the sums PostProcess(qroot_); // normalize densities @@ -787,27 +840,27 @@ class FastKde { // construct query and reference trees fx_timer_start(NULL, "tree_d"); - rroot_ = tree::MakeKdTreeMidpoint(rset_, leaflen, NULL); + rroot_ = tree::MakeKdTreeMidpoint(rset_, leaflen, + &r_old_from_new_); if(!strcmp(qfname, rfname)) { qset_.Alias(rset_); qroot_ = rroot_; + q_old_from_new_.Copy(r_old_from_new_); } else { Dataset query_dataset; query_dataset.InitFromFile(qfname); qset_.Own(&(query_dataset.matrix())); - qroot_ = tree::MakeKdTreeMidpoint(qset_, leaflen, NULL); + qroot_ = tree::MakeKdTreeMidpoint(qset_, leaflen, + &q_old_from_new_); } fx_timer_stop(NULL, "tree_d"); // initialize the density lists densities_l_.Init(qset_.n_cols()); - densities_l_.SetZero(); densities_e_.Init(qset_.n_cols()); - densities_e_.SetZero(); densities_u_.Init(qset_.n_cols()); - densities_u_.SetAll(rset_.n_cols()); // initialize the kernel kernel_.Init(fx_param_double_req(NULL, "bandwidth")); @@ -824,7 +877,9 @@ class FastKde { if((fname = fx_param_str(NULL, "fast_kde_output", NULL)) != NULL) { stream = fopen(fname, "w+"); } - densities_e_.PrintDebug("Fast KDE results", stream); + for(index_t q = 0; q < qset_.n_cols(); q++) { + fprintf(stream, "%g\n", densities_e_[q_old_from_new_[q]]); + } if(stream != stdout) { fclose(stream); diff --git a/fastlib/u/dongryel/kde/main.cc b/fastlib/u/dongryel/kde/main.cc index 2121253a32..785fe9ecb9 100644 --- a/fastlib/u/dongryel/kde/main.cc +++ b/fastlib/u/dongryel/kde/main.cc @@ -10,16 +10,25 @@ int main(int argc, char *argv[]) { FastKde fast_kde; fast_kde.Init(); fast_kde.Compute(fx_param_double(NULL, "tau", 0.1)); - //fast_kde.PrintDebug(); + + if(fx_param_exists(NULL, "fast_kde_output")) { + fast_kde.PrintDebug(); + } + Vector fast_kde_results; fast_kde_results.Alias(fast_kde.get_density_estimates()); if(do_naive) { NaiveKde naive_kde; - naive_kde.Init(fast_kde.get_query_dataset(), - fast_kde.get_reference_dataset()); + naive_kde.Init(fast_kde.get_query_dataset(), + fast_kde.get_query_old_from_new_mapping(), + fast_kde.get_reference_dataset(), + fast_kde.get_reference_old_from_new_mapping()); naive_kde.Compute(); - //naive_kde.PrintDebug(); + + if(fx_param_exists(NULL, "naive_kde_output")) { + naive_kde.PrintDebug(); + } naive_kde.ComputeMaximumRelativeError(fast_kde_results); }