From f5b7b3ee78e19e342d40130800ba8c82ffeedfbc Mon Sep 17 00:00:00 2001 From: rriegel Date: Thu, 30 Aug 2007 00:20:18 +0000 Subject: [PATCH] first stab at bichromatic/LOO --- fastlib/u/rriegel/nbc/nbc.cc | 347 ++++++++++++++++++++++++++++------- 1 file changed, 283 insertions(+), 64 deletions(-) diff --git a/fastlib/u/rriegel/nbc/nbc.cc b/fastlib/u/rriegel/nbc/nbc.cc index b4cad0c34d..77e8b43577 100644 --- a/fastlib/u/rriegel/nbc/nbc.cc +++ b/fastlib/u/rriegel/nbc/nbc.cc @@ -14,6 +14,7 @@ class Nbc { /** Point data includes class (referencs) and prior (queries). */ class NbcPoint { private: + index_t index_; Vector vec_; /** The point's class (if reference). */ bool is_pos_; @@ -26,6 +27,12 @@ class Nbc { public: + /** + * Gets the index. + */ + index_t index() const { + return index_; + } /** * Gets the vector. */ @@ -78,7 +85,8 @@ class Nbc { * @param data the vector data read from file */ template - void Set(const Param& param, const Vector& data) { + void Set(const Param& param, index_t index, const Vector& data) { + index_ = index; mem::Copy(vec_.ptr(), data.ptr(), vec_.length()); is_pos_ = (data[data.length() - 2] != 0.0); pi_ = data[data.length() - 1]; @@ -108,6 +116,9 @@ class Nbc { */ DRange const_pos; DRange const_neg; + /** Like const_pos but leaving out one point. Similar for _neg. */ + DRange const_pos_loo; + DRange const_neg_loo; /** The kernel for positive points. Similar for _neg. */ Kernel kernel_pos; Kernel kernel_neg; @@ -120,10 +131,17 @@ class Nbc { index_t count_neg; /** The specified threshold for certainty of positive class. */ double threshold; + /** Peak kernel values (to eliminate recomputation). */ + double peak_neg; + double peak_pos; + /** Whether to compute densities leaving out matching indices. */ + bool loo; OT_DEF_BASIC(Param) { OT_MY_OBJECT(const_pos); OT_MY_OBJECT(const_neg); + OT_MY_OBJECT(const_pos_loo); + OT_MY_OBJECT(const_neg_loo); OT_MY_OBJECT(kernel_pos); OT_MY_OBJECT(kernel_neg); OT_MY_OBJECT(dim); @@ -131,6 +149,9 @@ class Nbc { OT_MY_OBJECT(count_pos); OT_MY_OBJECT(count_neg); OT_MY_OBJECT(threshold); + OT_MY_OBJECT(peak_pos); + OT_MY_OBJECT(peak_neg); + OT_MY_OBJECT(loo); } public: @@ -141,11 +162,21 @@ class Nbc { kernel_pos.Init(fx_param_double_req(module, "h_pos")); kernel_neg.Init(fx_param_double_req(module, "h_neg")); threshold = fx_param_double(module, "threshold", 0.5); + loo = !fx_param_exists(module, "q"); } + /** + * Reflect some aspects of the data (Req THOR). + * + * Note: called *after* reading all data. + */ void SetDimensions(index_t vector_dimension, index_t n_points) { dim = vector_dimension; // last two cols already trimmed count_all = n_points; + peak_pos = kernel_pos.EvalUnnormOnSq(0) + / kernel_pos.CalcNormConstant(dim); + peak_neg = kernel_neg.EvalUnnormOnSq(0) + / kernel_neg.CalcNormConstant(dim); } /** @@ -154,15 +185,24 @@ class Nbc { void ComputeConsts(int count_pos_in, int count_neg_in) { double epsilon = min(threshold, 1 - threshold) * 1e-3; - double norm_pos = kernel_pos.CalcNormConstant(dim) * count_pos_in; + count_pos = count_pos_in; + count_neg = count_neg_in; + + double norm_pos = kernel_pos.CalcNormConstant(dim) * count_pos; const_pos.lo = (1 - threshold - epsilon) / norm_pos; const_pos.hi = (1 - threshold + epsilon) / norm_pos; - count_pos = count_pos_in; - double norm_neg = kernel_neg.CalcNormConstant(dim) * count_neg_in; + double norm_neg = kernel_neg.CalcNormConstant(dim) * count_neg; const_neg.lo = (threshold - epsilon) / norm_neg; const_neg.hi = (threshold + epsilon) / norm_neg; - count_neg = count_neg_in; + + norm_pos = kernel_pos.CalcNormConstant(dim) * (count_pos - 1); + const_pos_loo.lo = (1 - threshold - epsilon) / norm_pos; + const_pos_loo.hi = (1 - threshold + epsilon) / norm_pos; + + norm_neg = kernel_neg.CalcNormConstant(dim) * (count_neg - 1); + const_neg_loo.lo = (threshold - epsilon) / norm_neg; + const_neg_loo.hi = (threshold + epsilon) / norm_neg; ot::Print(dim); ot::Print(count_all); @@ -294,6 +334,8 @@ class Nbc { /** Bounds for query priors. Similar for _neg. */ DRange pi_pos; DRange pi_neg; + /** Bounds for query self-contrib (undone in LOO case). */ + DRange loo_contrib; OT_DEF_BASIC(NbcStat) { OT_MY_OBJECT(moment_info_pos); @@ -302,6 +344,9 @@ class Nbc { OT_MY_OBJECT(bound_neg); OT_MY_OBJECT(count_pos); OT_MY_OBJECT(count_neg); + OT_MY_OBJECT(pi_pos); + OT_MY_OBJECT(pi_neg); + OT_MY_OBJECT(loo_contrib); } public: @@ -319,6 +364,11 @@ class Nbc { count_neg = 0; pi_pos.InitEmptySet(); pi_neg.InitEmptySet(); + if (param.loo) { + loo_contrib.InitEmptySet(); + } else { + loo_contrib.Init(0, 0); + } } /** @@ -336,6 +386,13 @@ class Nbc { } pi_pos |= point.pi_pos(); pi_neg |= point.pi_neg(); + if (param.loo) { + if (point.is_pos()) { + loo_contrib |= param.peak_pos; + } else { + loo_contrib |= param.peak_neg; + } + } } /** @@ -351,6 +408,7 @@ class Nbc { count_neg += stat.count_neg; pi_pos |= stat.pi_pos; pi_neg |= stat.pi_neg; + loo_contrib |= stat.loo_contrib; } /** @@ -459,12 +517,33 @@ class Nbc { void Postprocess(const Param& param, const QPoint& q, index_t q_index, const RNode& r_root) { - if (param.const_pos.lo * density_pos * q.pi_pos() - > param.const_neg.hi * density_neg * q.pi_neg()) { - label &= LAB_POS; - } else if (param.const_neg.lo * density_neg * q.pi_neg() - > param.const_pos.hi * density_pos * q.pi_pos()) { - label &= LAB_NEG; + if (param.loo) { + // Withhold contribution of q from density for its class + if (q.is_pos()) { + if (param.const_pos_loo.lo * (density_pos - param.peak_pos) * q.pi_pos() + > param.const_neg.hi * density_neg * q.pi_neg()) { + label &= LAB_POS; + } else if (param.const_neg.lo * density_neg * q.pi_neg() + > param.const_pos_loo.hi * (density_pos - param.peak_pos) * q.pi_pos()) { + label &= LAB_NEG; + } + } else { + if (param.const_pos.lo * density_pos * q.pi_pos() + > param.const_neg_loo.hi * (density_neg - param.peak_neg) * q.pi_neg()) { + label &= LAB_POS; + } else if (param.const_neg_loo.lo * (density_neg - param.peak_neg) * q.pi_neg() + > param.const_pos.hi * density_pos * q.pi_pos()) { + label &= LAB_NEG; + } + } + } else { + if (param.const_pos.lo * density_pos * q.pi_pos() + > param.const_neg.hi * density_neg * q.pi_neg()) { + label &= LAB_POS; + } else if (param.const_neg.lo * density_neg * q.pi_neg() + > param.const_pos.hi * density_pos * q.pi_pos()) { + label &= LAB_NEG; + } } DEBUG_ASSERT_MSG(label != LAB_NEITHER, "Conflicting labels: [%g, %g]; %g > %g; %g > %g", @@ -766,14 +845,39 @@ class Nbc { DRange total_density_neg = unapplied_summary_results.density_neg + q_result->density_neg; - if (unlikely( - param.const_pos.lo * total_density_pos.lo * q.pi_pos() - > param.const_neg.hi * total_density_neg.hi * q.pi_neg())) { - q_result->label &= LAB_POS; - } else if (unlikely( - param.const_neg.lo * total_density_neg.lo * q.pi_neg() - > param.const_pos.hi * total_density_pos.hi * q.pi_pos())) { - q_result->label &= LAB_NEG; + if (param.loo) { + // Withhold contribution of q from density for its class + if (q.is_pos()) { + if (unlikely( + param.const_pos_loo.lo * (total_density_pos.lo - param.peak_pos) * q.pi_pos() + > param.const_neg.hi * total_density_neg.hi * q.pi_neg())) { + q_result->label &= LAB_POS; + } else if (unlikely( + param.const_neg.lo * total_density_neg.lo * q.pi_neg() + > param.const_pos_loo.hi * (total_density_pos.hi - param.peak_pos) * q.pi_pos())) { + q_result->label &= LAB_NEG; + } + } else { + if (unlikely( + param.const_pos.lo * total_density_pos.lo * q.pi_pos() + > param.const_neg_loo.hi * (total_density_neg.hi - param.peak_neg) * q.pi_neg())) { + q_result->label &= LAB_POS; + } else if (unlikely( + param.const_neg_loo.lo * (total_density_neg.lo - param.peak_neg) * q.pi_neg() + > param.const_pos.hi * total_density_pos.hi * q.pi_pos())) { + q_result->label &= LAB_NEG; + } + } + } else { + if (unlikely( + param.const_pos.lo * total_density_pos.lo * q.pi_pos() + > param.const_neg.hi * total_density_neg.hi * q.pi_neg())) { + q_result->label &= LAB_POS; + } else if (unlikely( + param.const_neg.lo * total_density_neg.lo * q.pi_neg() + > param.const_pos.hi * total_density_pos.hi * q.pi_pos())) { + q_result->label &= LAB_NEG; + } } } }; @@ -794,7 +898,7 @@ class Nbc { Delta* delta, GlobalResult* global_result, QPostponed* q_postponed) { - DEBUG_MSG(1.0, "tkde: ConsiderPairIntrinsic"); + DEBUG_MSG(1.0, "nbc: ConsiderPairIntrinsic"); double d_density_pos_hi = 0; if (r_node.stat().count_pos > 0) { @@ -813,7 +917,7 @@ class Nbc { } if (d_density_pos_hi == 0 && d_density_neg_hi == 0) { - DEBUG_MSG(1.0, "tkde: Exclusion"); + DEBUG_MSG(1.0, "nbc: Exclusion"); return false; } @@ -830,7 +934,7 @@ class Nbc { if (r_node.stat().count_neg > 0) { q_postponed->moment_info_neg.Add(r_node.stat().moment_info_neg); } - DEBUG_MSG(1.0, "tkde: Inclusion"); + DEBUG_MSG(1.0, "nbc: Inclusion"); return false; } #else @@ -842,13 +946,13 @@ class Nbc { if (r_node.stat().count_neg > 0) { q_postponed->moment_info_neg.Add(r_node.stat().moment_info_neg); } - DEBUG_MSG(1.0, "tkde: Inclusion"); + DEBUG_MSG(1.0, "nbc: Inclusion"); return false; } #endif - delta->d_density_pos.Init(0, r_node.stat().count_pos * d_density_pos_hi); - delta->d_density_pos.hi = + delta->d_density_pos.lo = 0; + delta->d_density_pos.hi = r_node.stat().count_pos * d_density_pos_hi; delta->d_density_neg.lo = 0; delta->d_density_neg.hi = r_node.stat().count_neg * d_density_neg_hi; @@ -872,31 +976,95 @@ class Nbc { const QSummaryResult& q_summary_result, const GlobalResult& global_result, QPostponed* q_postponed) { - DEBUG_ASSERT(q_summary_result.density_pos.lo < q_summary_result.density_pos.hi); - DEBUG_ASSERT(q_summary_result.density_neg.lo < q_summary_result.density_neg.hi); + DEBUG_ASSERT(q_summary_result.density_pos.lo <= q_summary_result.density_pos.hi); + DEBUG_ASSERT(q_summary_result.density_neg.lo <= q_summary_result.density_neg.hi); if (unlikely(q_summary_result.label != LAB_EITHER)) { DEBUG_ASSERT((q_summary_result.label & q_postponed->label) != LAB_NEITHER); q_postponed->label = q_summary_result.label; return false; - } else if (unlikely( - param.const_pos.lo - * q_summary_result.density_pos.lo - * q_node.stat().pi_pos.lo - > param.const_neg.hi - * q_summary_result.density_neg.hi - * q_node.stat().pi_neg.hi)) { - q_postponed->label = LAB_POS; - return false; - } else if (unlikely( - param.const_neg.lo - * q_summary_result.density_neg.lo - * q_node.stat().pi_neg.lo - > param.const_pos.hi - * q_summary_result.density_pos.hi - * q_node.stat().pi_pos.hi)) { - q_postponed->label = LAB_NEG; - return false; + } + + if (param.loo) { +#ifdef CHECK_POS_NEG_LOO_BOUNDS + // Make sure bounds hold for both pos and neg queries + if (unlikely( + (q_node.stat().count_pos == 0 || + param.const_pos_loo.lo + * (q_summary_result.density_pos.lo - param.peak_pos) + * q_node.stat().pi_pos.lo + > param.const_neg.hi + * q_summary_result.density_neg.hi + * q_node.stat().pi_neg.hi) + && (q_node.stat().count_neg == 0 || + param.const_pos.lo + * q_summary_result.density_pos.lo + * q_node.stat().pi_pos.lo + > param.const_neg_loo.hi + * (q_summary_result.density_neg.hi - param.peak_neg) + * q_node.stat().pi_neg.hi))) { + q_postponed->label = LAB_POS; + return false; + } else if (unlikely( + (q_node.stat().count_pos == 0 || + param.const_neg.lo + * q_summary_result.density_neg.lo + * q_node.stat().pi_neg.lo + > param.const_pos_loo.hi + * (q_summary_result.density_pos.hi - param.peak_pos) + * q_node.stat().pi_pos.hi) + && (q_node.stat().count_neg == 0 || + param.const_neg_loo.lo + * (q_summary_result.density_neg.lo - param.peak_neg) + * q_node.stat().pi_neg.lo + > param.const_pos.hi + * q_summary_result.density_pos.hi + * q_node.stat().pi_pos.hi))) { + q_postponed->label = LAB_NEG; + return false; + } +#else + // Note const_pos.lo < const_pos_loo.lo, etc. + if (unlikely( + param.const_pos.lo + * (q_summary_result.density_pos.lo - param.peak_pos) + * q_node.stat().pi_pos.lo + > param.const_neg_loo.hi + * q_summary_result.density_neg.hi + * q_node.stat().pi_neg.hi)) { + q_postponed->label = LAB_POS; + return false; + } else if (unlikely( + param.const_neg.lo + * (q_summary_result.density_neg.lo - param.peak_neg) + * q_node.stat().pi_neg.lo + > param.const_pos_loo.hi + * q_summary_result.density_pos.hi + * q_node.stat().pi_pos.hi)) { + q_postponed->label = LAB_NEG; + return false; + } +#endif + } else { + if (unlikely( + param.const_pos.lo + * q_summary_result.density_pos.lo + * q_node.stat().pi_pos.lo + > param.const_neg.hi + * q_summary_result.density_neg.hi + * q_node.stat().pi_neg.hi)) { + q_postponed->label = LAB_POS; + return false; + } else if (unlikely( + param.const_neg.lo + * q_summary_result.density_neg.lo + * q_node.stat().pi_neg.lo + > param.const_pos.hi + * q_summary_result.density_pos.hi + * q_node.stat().pi_pos.hi)) { + q_postponed->label = LAB_NEG; + return false; + } } return true; @@ -917,16 +1085,28 @@ class Nbc { }; void NbcMain(datanode *module) { + // TODO: LOO, reporting, multi-bw, recursive "bfs", multi-thresh + // make sure bichromatic works + // for LOO: easiest to include self but correct when pruning + // i.e. subtract lb self-contrib from ub dens; vice versa + // make sure to use correct coefficients + // lb dens - ub self-contrib should be clamped positive + // lb/ub self-contrib deps on bws of pos/neg, what's in node + //thor::MonochromaticDualTreeMain >( // module, "tkde"); - const char *gnp_name = "tkde"; + const char *gnp_name = "nbc"; const int DATA_CHANNEL = 110; const int Q_RESULTS_CHANNEL = 120; const int GNP_CHANNEL = 200; double results_megs = fx_param_double(module, "results/megs", 1000); - DistributedCache *points_cache; - index_t n_points; - ThorTree tree; DistributedCache q_results; + DistributedCache *q_points_cache; + DistributedCache *r_points_cache; + index_t n_q_points; + index_t n_r_points; + ThorTree *q_tree; + ThorTree *r_tree; + DistributedCache q_results; Nbc::Param param; rpc::Init(); @@ -940,38 +1120,77 @@ void NbcMain(datanode *module) { param.Init(fx_submodule(module, gnp_name, gnp_name)); fx_timer_start(module, "read"); - points_cache = new DistributedCache(); - n_points = thor::ReadPoints( + r_points_cache = new DistributedCache(); + n_r_points = thor::ReadPoints( param, DATA_CHANNEL + 0, DATA_CHANNEL + 1, - fx_submodule(module, "data", "data"), points_cache); + fx_submodule(module, "r", "r"), + r_points_cache); + if (fx_param_exists(module, "q")) { + q_points_cache = new DistributedCache(); + n_q_points = thor::ReadPoints( + param, DATA_CHANNEL + 2, DATA_CHANNEL + 3, + fx_submodule(module, "q", "q"), + q_points_cache); + } else { + q_points_cache = r_points_cache; + n_q_points = n_r_points; + } fx_timer_stop(module, "read"); - Nbc::QPoint default_point; - CacheArray::GetDefaultElement( - points_cache, &default_point); - param.SetDimensions(default_point.vec().length(), n_points); + Nbc::RPoint default_point; + CacheArray::GetDefaultElement( + r_points_cache, &default_point); + param.SetDimensions(default_point.vec().length(), n_r_points); fx_timer_start(module, "tree"); - thor::CreateKdTree( - param, DATA_CHANNEL + 2, DATA_CHANNEL + 3, - fx_submodule(module, "tree", "tree"), n_points, points_cache, &tree); + r_tree = new ThorTree(); + thor::CreateKdTree( + param, DATA_CHANNEL + 4, DATA_CHANNEL + 5, + fx_submodule(module, "r_tree", "r_tree"), + n_r_points, r_points_cache, r_tree); + if (fx_param_exists(module, "q")) { + q_tree = new ThorTree(); + thor::CreateKdTree( + param, DATA_CHANNEL + 6, DATA_CHANNEL + 7, + fx_submodule(module, "q_tree", "q_tree"), + n_q_points, q_points_cache, q_tree); + } else { + q_tree = r_tree; + } fx_timer_stop(module, "tree"); // This should have been a first-order reduce at the time of read - param.ComputeConsts(tree.root().stat().count_pos, - tree.root().stat().count_neg); + param.ComputeConsts(r_tree->root().stat().count_pos, + r_tree->root().stat().count_neg); Nbc::QResult default_result; default_result.Init(param); - tree.CreateResultCache(Q_RESULTS_CHANNEL, default_result, - results_megs, &q_results); + q_tree->CreateResultCache(Q_RESULTS_CHANNEL, default_result, + results_megs, &q_results); Nbc::GlobalResult *global_result; thor::RpcDualTree >( fx_submodule(module, "gnp", "gnp"), GNP_CHANNEL, param, - &tree, &tree, &q_results, &global_result); + q_tree, r_tree, &q_results, &global_result); delete global_result; + // Emit the results; this needs to be folded into THOR + Matrix classifications; + classifications.Init(1, n_q_points); + if (rpc::is_root()) { + CacheArray result_array; + CacheArray points_array; + result_array.Init(&q_results, BlockDevice::M_READ); + points_array.Init(q_points_cache, BlockDevice::M_READ); + CacheReadIter result_iter(&result_array, 0); + CacheReadIter points_iter(&points_array, 0); + for (index_t i = 0; i < n_q_points; i++, + result_iter.Next(), points_iter.Next()) { + classifications.set(0, (*points_iter).index(), (*result_iter).label); + } + } + data::Save(fx_param_str(module, "out", "out.csv"), classifications); + rpc::Done(); }