apaam/mlpack
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

first stab at bichromatic/LOO

Browse Source
This commit is contained in:
rriegel
2007-08-30 00:20:18 +00:00
parent 75cb14ae1f
commit f5b7b3ee78
1 changed files with 283 additions and 64 deletions
Download Patch File Download Diff File Expand all files Collapse all files
fastlib/u/rriegel/nbc/nbc.cc
+283 -64
View File
@@ -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<typename Param>
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<Tkde, DualTreeDepthFirst<Tkde> >(
// 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<Nbc::Param, Nbc::QPoint, Nbc::QNode> tree; DistributedCache q_results;
DistributedCache *q_points_cache;
DistributedCache *r_points_cache;
index_t n_q_points;
index_t n_r_points;
ThorTree<Nbc::Param, Nbc::QPoint, Nbc::QNode> *q_tree;
ThorTree<Nbc::Param, Nbc::RPoint, Nbc::RNode> *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<Nbc::QPoint>(
r_points_cache = new DistributedCache();
n_r_points = thor::ReadPoints<Nbc::RPoint>(
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<Nbc::QPoint>(
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<Nbc::QPoint>::GetDefaultElement(
points_cache, &default_point);
param.SetDimensions(default_point.vec().length(), n_points);
Nbc::RPoint default_point;
CacheArray<Nbc::RPoint>::GetDefaultElement(
r_points_cache, &default_point);
param.SetDimensions(default_point.vec().length(), n_r_points);
fx_timer_start(module, "tree");
thor::CreateKdTree<Nbc::QPoint, Nbc::QNode>(
param, DATA_CHANNEL + 2, DATA_CHANNEL + 3,
fx_submodule(module, "tree", "tree"), n_points, points_cache, &tree);
r_tree = new ThorTree<Nbc::Param, Nbc::RPoint, Nbc::RNode>();
thor::CreateKdTree<Nbc::RPoint, Nbc::RNode>(
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<Nbc::Param, Nbc::QPoint, Nbc::QNode>();
thor::CreateKdTree<Nbc::QPoint, Nbc::QNode>(
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<Nbc, DualTreeDepthFirst<Nbc> >(
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<Nbc::QResult> result_array;
CacheArray<Nbc::QPoint> points_array;
result_array.Init(&q_results, BlockDevice::M_READ);
points_array.Init(q_points_cache, BlockDevice::M_READ);
CacheReadIter<Nbc::QResult> result_iter(&result_array, 0);
CacheReadIter<Nbc::QPoint> 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();
}