Algorithm-Classifier-IsolationForest
view release on metacpan or search on metacpan
lib/Algorithm/Classifier/IsolationForest.pm view on Meta::CPAN
use constant EULER => 0.5772156649015329;
# Narrowed to C double precision so _randn() multiplies by the exact
# constant _c_randn() uses. A no-op on nvsize == 8 perls.
use constant TWO_PI => unpack( 'd', pack 'd', 6.283185307179586 );
# Node-type tags stored in index 0 of every tree node arrayref.
# 0 is falsy, so while ($node->[0]) acts as while (!leaf).
use constant _NODE_LEAF => 0;
use constant _NODE_AXIS => 1;
use constant _NODE_OBLIQUE => 2;
# The Inline::C tree builder computes everything in C doubles. On a perl
# whose NV is wider than a double (-Duselongdouble / -Dusequadmath) the
# pure-Perl builder keeps extra low bits at every step, so the two
# backends would stop producing bit-identical trees for the same seed
# (the parity t/03-fit-determinism.t checks). _NV_IS_DOUBLE guards
# narrowing statements wherever the pure-Perl builder computes a value
# that gets STORED in a tree (split points, hyperplane coefficients and
# offsets, impute fills), rounding at the same points the C builder
# rounds. It is compile-time true on nvsize == 8 perls, so there the
# guarded statements are optimised away and cost nothing.
#
# The row-partition loops (v < split, dot <= b) are deliberately NOT
# narrowed: with both operands already double-exact an axis comparison
# is identical anyway, and an oblique dot product accumulated in a wider
# NV flips a comparison only when |dot - b| falls inside the NV-vs-double
# rounding gap (~1e-19 relative) -- negligible, and those are the hot
# loops.
use constant _NV_IS_DOUBLE => $Config::Config{nvsize} == 8;
# Round an NV to C double precision. Only ever reached on wide-NV
# perls -- see _NV_IS_DOUBLE.
sub _to_double { unpack 'd', pack 'd', $_[0] }
# ---------------------------------------------------------------------------
# Optional Inline::C accelerator for the scoring hot path.
#
# pack_input_xs(data_sv, out_sv, n_pts, n_feats, miss_mode, fill_sv)
# Walks the Perl arrayref-of-arrayrefs and writes a packed double buffer
# into out_sv. Replaces the dominant per-call Perl map-pack loop.
# miss_mode selects how an undef cell is packed: 0 => 0.0, 1 => the
# per-feature fill from fill_sv (impute), 2 => NaN (nan strategy).
#
# score_all_xs(nodes_av, idx_av, val_av, x_sv, sm_sv,
# n_pts, n_feats, n_trees, use_openmp)
# Sums path lengths for all n_pts query points across all n_trees trees
# in one call. Outer loop over points is OpenMP-parallel when the
# module was built with OpenMP (each iteration writes to a unique sm[i],
# so no synchronisation is needed). Tree pointers are extracted from
# the AVs before the parallel region; the parallel region touches only
# raw int / double buffers.
#
# vote_all_xs(nodes_av, idx_av, val_av, x_sv, sm_sv,
# n_pts, n_feats, n_trees, depth_cut, min_votes, use_openmp)
# Majority-voting (voting => 'majority') counterpart of score_all_xs:
# instead of summing path lengths it counts, per point, how many trees
# "vote anomalous" (path length <= depth_cut). min_votes == 0 writes
# the full vote count into sm[i]; min_votes > 0 writes a 0.0/1.0 label
# with per-point early exit once the majority outcome is decided --
# the MVIForest scoring loop. See the function's own comment.
#
# Node layout (6 doubles per node, "IF_NZ = 6"):
# leaf: [0, size, c(size), 0, 0, 0]
# axis: [1, attr, split, li, ri, 0]
# oblique: [2, coff, nf, li, ri, b]
#
# c(size) is the expected-path-length adjustment for a leaf holding
# `size` points, precomputed by _pack_tree (it involves a log(); doing
# it at pack time keeps transcendentals out of the per-point per-tree
# scoring loop). The fit-time TreeBuf writer leaves that slot 0 --
# its buffers are unpacked into Perl trees and re-packed by
# _pack_tree before score_all_xs ever sees them.
#
# Coefficient storage uses a Structure-of-Arrays layout: one int32 array
# per tree (feature indices, packed with 'l*') and one double array per
# tree (coefficients, packed with 'd*'). Both are indexed by `coff` --
# the same offset addresses paired entries in the two arrays. Splitting
# them this way halves index bandwidth, removes the per-element
# (int)<double> cast inside the SIMD loop, and lets the value loads be
# contiguous so the compiler emits a clean FMA chain over val[k] with
# the feature gather on xi[idx[k]] kept separate.
#
# Dense-pack fast path: when an oblique node uses every feature (the
# common case in extended mode with extension_level == n_features - 1),
# _pack_tree writes its coefficients in feature order so val[k] is the
# coefficient for feature k. score_all_xs detects this via `nf ==
# n_feats` and uses a no-gather dot product (dot += val[k] * xi[k])
# that vectorizes cleanly with FMA -- substantially faster than the
# sparse gather path on high-feature-count models.
# x: row-major doubles, n_pts rows of n_feats each.
# sums: out double array of length n_pts; score_all_xs writes once per i.
#
# OpenMP is enabled at module load when the toolchain accepts -fopenmp and
# libgomp is linkable; otherwise the same C code compiles to a serial loop
# (the #pragma is silently ignored without _OPENMP defined).
# ---------------------------------------------------------------------------
our $HAS_C = 0;
our $HAS_OPENMP = 0;
our $HAS_SIMD = 0;
our $OPT_LEVEL = ''; # the actual -O.../-march=... flags used to build, if any
our $C_SOURCE = ''; # 'prebuilt' (object installed at `make` time) or
# 'runtime' (compiled at first load into _Inline/);
# '' when $HAS_C is 0
{
my $C_CODE = <<'__INLINE_C__';
#include <math.h>
#include <string.h>
#include <stdint.h>
#ifdef _OPENMP
#include <omp.h>
#endif
#define IF_NZ 6
/* Data prefetch hint; a no-op on compilers without __builtin_prefetch.
* Purely a performance hint -- never affects results. */
#if defined(__GNUC__) || defined(__clang__)
#define IF_PREFETCH(p) __builtin_prefetch(p)
#else
#define IF_PREFETCH(p)
lib/Algorithm/Classifier/IsolationForest.pm view on Meta::CPAN
}
} else {
/* Sparse oblique split: only nf < n_feats features
* participate, so we still need the gather on
* xi[idx_p[k]]. Storing idx as contiguous int32
* (rather than interleaved doubles) keeps the gather
* pattern clean and the val[] load contiguous. */
const int *idx_p = ico + (size_t)coff;
#ifdef _OPENMP
#pragma omp simd reduction(+:dot)
#endif
for (int k = 0; k < nf; k++) {
dot += val_p[k] * xi[idx_p[k]];
}
}
ni = (dot <= b) ? li : ri;
}
depth++;
}
}
/* score_all_xs(nodes_av, idx_av, val_av, x_sv, sm_sv,
* n_pts, n_feats, n_trees, use_openmp)
*
* Scores all points across all trees in one C call. See header comment
* above for the bigger picture. Writes sm[i] = sum_over_trees(path_len);
* the caller need not zero-init sm.
*
* idx_av holds per-tree packed int32 buffers of feature indices and
* val_av holds per-tree packed double buffers of coefficients (the SoA
* counterpart of the old interleaved layout). See the file-top
* comment for the rationale.
*
* Thread-safety: the parallel region only reads node/idx/val/x pointers
* (extracted before the region) and writes sm[i] for a unique i per
* iteration. No Perl API is called from inside the parallel region. */
void score_all_xs(SV* nodes_av_sv, SV* idx_av_sv, SV* val_av_sv,
SV* x_sv, SV* sm_sv,
int n_pts, int n_feats, int n_trees,
int use_openmp){
STRLEN tl;
AV *nodes_av, *idx_av, *val_av;
const double *xd;
double *sm;
int ti;
if (!SvROK(nodes_av_sv) || SvTYPE(SvRV(nodes_av_sv)) != SVt_PVAV ||
!SvROK(idx_av_sv) || SvTYPE(SvRV(idx_av_sv)) != SVt_PVAV ||
!SvROK(val_av_sv) || SvTYPE(SvRV(val_av_sv)) != SVt_PVAV) {
croak("score_all_xs: nodes/idx/val must be arrayrefs");
}
nodes_av = (AV*)SvRV(nodes_av_sv);
idx_av = (AV*)SvRV(idx_av_sv);
val_av = (AV*)SvRV(val_av_sv);
/* C99 VLAs -- n_trees is small (typ. 100) and fits on the stack. */
const double *node_ptrs[n_trees];
const int *idx_ptrs[n_trees];
const double *val_ptrs[n_trees];
/* forest_bytes totals every buffer the tree walks touch; it decides
* between the two loop shapes below. */
size_t forest_bytes = 0;
for (ti = 0; ti < n_trees; ti++) {
SV** np = av_fetch(nodes_av, ti, 0);
SV** ip = av_fetch(idx_av, ti, 0);
SV** vp = av_fetch(val_av, ti, 0);
if (!np || !*np || !ip || !*ip || !vp || !*vp) {
croak("score_all_xs: missing tree %d", ti);
}
node_ptrs[ti] = (const double*)SvPVbyte(*np, tl); forest_bytes += tl;
idx_ptrs[ti] = (const int*) SvPVbyte(*ip, tl); forest_bytes += tl;
val_ptrs[ti] = (const double*)SvPVbyte(*vp, tl); forest_bytes += tl;
}
xd = (const double*)SvPVbyte(x_sv, tl);
sm = (double*)SvPVbyte_force(sm_sv, tl);
/* Two loop shapes over the same per-point ascending-t additions --
* bit-identical results either way, so the size heuristic choosing
* between them can never change scores.
*
* Point-major (small forests): each point walks all trees with its
* path-length sum held in a register. Cheapest per walk, and the
* whole forest stays cache-resident across points anyway.
*
* Tree-blocked (large forests): once the forest outgrows L3, the
* point-major loop re-streams every tree's nodes and coefficients
* from memory for every point -- an extended-mode tree is ~56 KB
* at 16 features (24 KB nodes + 32 KB dense coefficients), and its
* per-tree scoring cost measured 2.2x worse at 400 trees than at
* 100. Walking a block of points through ONE tree at a time keeps
* that tree hot in L1/L2 while the block's rows stream through it
* (measured 3.1x faster at 400 extended trees, 20k points). The
* blocked shape pays an sm[i] load+store per walk instead of a
* register add, which measurably hurts cheap axis walks while the
* forest still fits in cache -- hence the byte threshold rather
* than always tiling. */
if (forest_bytes <= (size_t)4 * 1024 * 1024) {
#ifdef _OPENMP
#pragma omp parallel for schedule(static) if(use_openmp)
#endif
for (int i = 0; i < n_pts; i++) {
const double *xi = xd + (size_t)i * (size_t)n_feats;
double sum = 0.0;
for (int t = 0; t < n_trees; t++) {
sum += if_walk_tree(node_ptrs[t], idx_ptrs[t],
val_ptrs[t], xi, n_feats);
}
sm[i] = sum;
}
}
else {
/* 256 rows x 16 features x 8 bytes = 32 KB of input per block
* -- comfortable in L2 next to one tree. Each OpenMP thread
* owns whole blocks and therefore a unique slice of sm[], so
* there is still no synchronisation. For small batches the
* tile shrinks to keep ~4 blocks per thread available; losing
* per-block tree reuse there is fine, since a small batch
* never re-streams much anyway. */
int tile = 256;
lib/Algorithm/Classifier/IsolationForest.pm view on Meta::CPAN
# throttling under sustained heavy use and can make throughput
# *worse* than a conservative target like x86-64-v3 (AVX2, no
# AVX-512). Either way, the cached artefact under _Inline/ is then
# pinned to that instruction set, so leave both unset if the
# directory is shared across machines with different CPUs.
my $arch = $def_arch;
if ( defined $ENV{IF_ARCH} ) {
if ( $ENV{IF_ARCH} eq '' or $ENV{IF_ARCH} eq 'none' ) {
# Explicit opt-out: overrides an arch recorded at configure
# time (there is no other way to request a plain build on
# an install configured with IF_ARCH).
$arch = '';
} elsif ( $ENV{IF_ARCH} =~ /\A[A-Za-z0-9_.+=-]+\z/ ) {
$arch = $ENV{IF_ARCH};
} else {
warn "Algorithm::Classifier::IsolationForest: ignoring invalid " . "IF_ARCH value '$ENV{IF_ARCH}'\n";
}
} elsif ( $ENV{IF_NATIVE} ) {
$arch = 'native';
}
# -ffp-contract=off rides along with any -march: once the target
# has FMA (x86-64-v3, most -march=native hosts), the compiler may
# otherwise contract a*b+c expressions into fused multiply-adds
# whose different rounding breaks the documented guarantee that
# use_c => 1 and use_c => 0 build bit-identical trees (one ulp in a
# split value cascades into a structurally different tree). The
# -march speedup comes from AVX2 vectorization, not contraction,
# so this costs little (verified against the fit-determinism and
# scoring-parity tests).
my $opt_level = $opt;
$opt_level .= " -march=$arch -ffp-contract=off" if length $arch;
# IF_NO_OPENMP=1 forces the serial C build: the OpenMP compile attempt
# is skipped, so the object has no libgomp linkage and never starts an
# OpenMP runtime in the process. Distinct from OMP_NUM_THREADS=1,
# which runs the parallel code on a single thread but still loads
# libgomp. An explicit IF_NO_OPENMP=0 re-enables OpenMP over a
# no-openmp configure-time default.
my $no_omp
= defined $ENV{IF_NO_OPENMP}
? ( $ENV{IF_NO_OPENMP} ? 1 : 0 )
: $def_no_omp;
# The prebuilt object is only trusted when the effective flags match
# what it was compiled with; any difference -- or an explicit
# IF_RUNTIME_BUILD=1 -- falls through to the classic runtime Inline::C
# build below, which honours the requested flags via the MD5-keyed
# _Inline/ cache exactly as before prebuilt support existed.
# IF_INSTALL_BUILD is the `make` rule driving the install-time compile
# (see Makefile.PL); it must never short-circuit into loading an
# older object.
my $use_prebuilt
= $prebuilt
&& !$ENV{IF_RUNTIME_BUILD}
&& !$ENV{IF_INSTALL_BUILD}
&& $opt eq $def_opt
&& $arch eq $def_arch
&& $no_omp == $def_no_omp;
# Inline::C hashes the C source to decide whether to rebuild but
# does NOT include CCFLAGS / OPTIMIZE in that hash. Without the
# tag below, toggling IF_NATIVE/IF_ARCH/IF_OPT (or editing the
# optimisation flags here) would silently reuse a cached binary
# built with stale flags. Embedding the active flags as a leading
# comment forces the hash to differ when they change. The OpenMP
# and serial builds get distinct tags so they cache to separate
# artefacts.
my $omp_tag = "/* if_build: openmp $opt_level */\n";
my $serial_tag = "/* if_build: serial $opt_level */\n";
if ( $ENV{IF_INSTALL_BUILD} ) {
# `make` is driving: the rule Makefile.PL appended runs this load
# with IF_INSTALL_BUILD=1 and @ARGV = (version, INST_ARCHLIB),
# which is where Inline's install mode reads them from. _INSTALL_
# makes Inline compile the backend and place the shared object
# under blib/arch so `make install` ships it; NAME/VERSION give
# the object a fixed identity XSLoader can find at run time
# (Inline's install mode also requires both and checks VERSION
# against $ARGV[0]). Same OpenMP-then-serial fallback as the
# runtime build below.
my @install = (
NAME => __PACKAGE__,
VERSION => $VERSION,
_INSTALL_ => 1,
);
unless ($no_omp) {
local $@;
eval {
require Inline;
Inline->import(
C => $omp_tag . $C_CODE,
CCFLAGS => '-fopenmp',
OPTIMIZE => $opt_level,
LIBS => '-lm -lgomp',
@install,
);
$HAS_C = 1;
};
} ## end unless ($no_omp)
unless ($HAS_C) {
local $@;
eval {
require Inline;
Inline->import(
C => $serial_tag . $C_CODE,
OPTIMIZE => $opt_level,
LIBS => '-lm',
@install,
);
$HAS_C = 1;
};
} ## end unless ($HAS_C)
$C_SOURCE = 'prebuilt' if $HAS_C;
} else {
# Fast path: the object compiled at `make` time was installed
# under auto/ like any XS module, so plain XSLoader digs it out of
# @INC with no Inline involvement -- no compiler, no _Inline/
# directory, and a few ms instead of a first-run compile. Any
lib/Algorithm/Classifier/IsolationForest.pm view on Meta::CPAN
default :: undef
- seed :: optional integer to seed srand with for reproducible trees...
see perldoc -f srand for more info. This number is processed via abs(int()).
default :: undef
- mode :: if it should be IF or EIF
axis :: classic axis-parallel splits (IF)
extended :: oblique hyperplane splits (EIF)
default :: axis
- extension_level :: extended mode only... how many features take partin each
split. 0 behaves like a single-feature (axis) cut; the
maximum (n_features - 1) uses every varying feature. undef
=> maximum. Clamped to [0, n_features - 1] at fit time.
- contamination :: expected fraction of anomalies, in (0, 0.5]. When given,
fit() learns a score threshold that flags this fraction of
the training set, and predict() uses it by default. undef
=> no learned threshold (predict() falls back to 0.5).
default :: undef
- missing :: how fit() treats undef (missing) feature cells. Scoring always
tolerates undef regardless of this setting; it governs fit().
die :: croak from fit() if the training data contains any
undef cell. Scoring still maps undef to 0 (the
long-standing behaviour), so a model fitted on clean
data can still score rows with missing features.
zero :: treat a missing cell as the value 0, at fit and score.
impute :: replace a missing cell with the per-feature mean (or
median, see impute_with) learned from the present
values at fit time. The fill vector is stored on the
model and reused for scoring and persistence.
nan :: build feature ranges from present values only and route
a point missing the split feature to the right child,
consistently at fit and score time. Missingness is
preserved as signal rather than filled.
default :: die
- impute_with :: 'mean' or 'median'; the statistic used to compute the
per-feature fill under missing => 'impute'. Ignored otherwise.
default :: mean
- voting :: how the per-tree results are aggregated at scoring time.
Trees are built identically in both settings -- only aggregation
changes -- so the knob composes with either mode (axis or
extended) and an existing model may switch it after the fact with
set_voting() (which relearns a contamination threshold for the
new mode).
mean :: classic Isolation Forest: a sample's path lengths
across all trees are averaged and normalised into
one anomaly score; predict() thresholds that score.
majority :: Majority Voting Isolation Forest (MVIForest;
Chabchoub, Togbe, Boly & Chiky 2022 -- see
REFERENCES). Each tree scores the sample on its own
(s_i = 2**(-h_i / c(psi))) and votes it anomalous
when s_i >= the decision threshold; predict() flags
the sample when more than half of the trees
(int(n_trees/2) + 1) vote anomalous, and stops
walking trees per sample as soon as the outcome is
decided. The threshold argument/default of the
predict methods is therefore the PER-TREE cutoff
here, not a forest-level score cutoff.
score_samples() returns the fraction of trees
voting anomalous -- still in [0, 1], but discrete
in steps of 1/n_trees. contamination composes: fit()
learns the per-tree cutoff that flags the requested
fraction of the training set.
default :: mean
- parallel_fit :: positive integer N => build the trees across N forked
worker processes during fit(). Each worker gets a derived seed
(parent seed + worker_id * 1009) so the parallel fit is
reproducible across runs at fixed worker count -- but the trees
produced are NOT bit-identical to a serial fit with the same
seed, because the RNG draws happen in a different order.
Inference is unaffected. Falls back silently to serial on
platforms without a real fork() (e.g. Windows without Cygwin).
default :: undef (serial)
- use_c :: boolean, override whether the Inline::C backend is used for
both scoring and fit()'s tree builder. When false the instance
falls back to pure Perl for both even if the C backend compiled
successfully. When true (or unset) the C backend is used if
available ($HAS_C). fit() with use_c on produces bit-identical
trees to use_c off for the same seed -- only build speed differs.
default :: $HAS_C
- use_openmp :: boolean, override whether OpenMP parallel scoring is
used inside score_all_xs(). When false the C tree walk runs
single-threaded even if OpenMP was linked in. Ignored when
use_c is false (pure Perl has no OpenMP path).
default :: $HAS_OPENMP
- use_openmp_fit :: boolean, build fit()'s trees across OpenMP threads
(one tree per thread) instead of the single-threaded C builder.
Opt-in and off by default: unlike use_c/use_openmp, this changes
which trees get built. Perl's RNG isn't safe to call from
multiple OS threads sharing one interpreter, so this path seeds
an independent PRNG per tree from the tree index rather than
Drand01() -- trees differ from the use_c (single-threaded)
and pure-Perl paths even with the same seed, though a fixed
seed and n_trees still reproduce the same trees regardless of
OMP_NUM_THREADS or scheduling. Does NOT compose with
parallel_fit: a forked child starting its own OpenMP region
after the parent process has used OpenMP for anything can
hang (a general fork()+libgomp limitation), so parallel_fit's
workers always use the single-threaded C builder regardless
of this setting -- setting both just means parallel_fit wins.
Ignored (clamped to 0) when use_c is false or OpenMP isn't
linked in.
default :: 0
Note: log2 under Perl is as below...
log($psi) / log(2)
=cut
sub new {
my ( $class, %args ) = @_;
lib/Algorithm/Classifier/IsolationForest.pm view on Meta::CPAN
=head2 path_lengths(\@data)
Returns an arrayref of the mean isolation depth per sample, for inspection.
my $lengths = $forest->path_lengths(\@data);
print "x, y, length\n";
my $int=0;
while (defined($data[$int])) {
print $data[$int][0].', '.$data[$int][1].', '.$lengths->[$int]."\n";
$int++;
}
=cut
sub path_lengths {
my ( $self, $data ) = @_;
$self->_check_fitted;
my $trees = $self->{trees};
my $t = scalar @$trees;
if ( $self->{_use_c} && $self->{_c_nodes} ) {
my ( $n_pts, $nf, $x_packed ) = $self->_resolve_input($data);
my $sums_packed = "\0" x ( $n_pts * 8 );
score_all_xs(
$self->{_c_nodes}, $self->{_c_coef_idx}, $self->{_c_coef_val},
$x_packed, $sums_packed, $n_pts,
$nf, $t, $self->{_use_openmp}
);
my $result = [];
finalize_path_lengths_xs( $sums_packed, $n_pts, $t + 0.0, $result );
return $result;
} ## end if ( $self->{_use_c} && $self->{_c_nodes} )
$data = $self->_prepare_perl_input($data);
my $nan = $self->{missing} eq 'nan' ? 1 : 0;
# Pure-Perl fallback (tree-outer, sample-inner for cache locality).
my @sums = (0) x @$data;
for my $tree (@$trees) {
for my $i ( 0 .. $#$data ) {
$sums[$i] += _path_length( $data->[$i], $tree, 0, $nan );
}
}
return [ map { $_ / $t } @sums ];
} ## end sub path_lengths
=head2 predict(\@data, $threshold)
Returns an arrayref of 0/1 labels for the specified data.
If threshold is not specified it uses the contamination-learned cutoff (if
C<fit> was called with C<contamination>), otherwise 0.5.
Under C<< voting => 'majority' >> the threshold is the per-tree score
cutoff each tree votes against, and a sample is labelled 1 when more than
half of the trees (C<int(n_trees/2) + 1>) vote it anomalous. Tree walking
stops per sample as soon as the outcome is decided, so this is typically
cheaper than scoring.
my $results = $forest->predict(\@data, $threshold);
print "x, y, result\n";
my $int=0;
while (defined($data[$int])) {
print $data[$int][0].', '.$data[$int][1].', '.$results->[$int]."\n";
$int++;
}
=cut
sub predict {
my ( $self, $data, $threshold ) = @_;
$threshold
= defined $threshold ? $threshold
: defined $self->{threshold} ? $self->{threshold}
: 0.5;
$self->_check_fitted;
# Majority voting: $threshold is the PER-TREE score cutoff and the
# label is the majority of the tree votes (int(t/2) + 1). Both the C
# and the Perl loop stop walking a sample's remaining trees as soon
# as the outcome is decided -- MVIForest's "stop at majority" saving.
if ( $self->{voting} eq 'majority' ) {
my $trees = $self->{trees};
my $t = scalar @$trees;
my $cut = _depth_cut( $threshold, $self->{c_psi} );
my $maj = _min_votes($t);
if ( $self->{_use_c} && $self->{_c_nodes} ) {
my ( $n_pts, $nf, $x_packed ) = $self->_resolve_input($data);
my $labels_packed = "\0" x ( $n_pts * 8 );
vote_all_xs( $self->{_c_nodes}, $self->{_c_coef_idx}, $self->{_c_coef_val},
$x_packed, $labels_packed, $n_pts, $nf, $t, $cut, $maj, $self->{_use_openmp} );
my $result = [];
vote_labels_xs( $labels_packed, $n_pts, $result );
return $result;
}
my $rows = $self->_prepare_perl_input($data);
my $nan = $self->{missing} eq 'nan' ? 1 : 0;
my @labels;
for my $x (@$rows) {
my $votes = 0;
my $label = 0;
for my $ti ( 0 .. $t - 1 ) {
if ( _path_length( $x, $trees->[$ti], 0, $nan ) <= $cut ) {
$votes++;
if ( $votes >= $maj ) { $label = 1; last }
}
last if $votes + ( $t - 1 - $ti ) < $maj;
}
push @labels, $label;
} ## end for my $x (@$rows)
return \@labels;
} ## end if ( $self->{voting} eq 'majority' )
# Fast path: threshold the raw path-length sums directly, skipping the
# per-point exp() and the intermediate scores arrayref.
# Derivation: score = exp(-sum * log(2) / (c*t))
# so score >= T iff sum <= -log(T) * c * t / log(2)
# Only valid for a normal threshold in (0, 1) and a positive c.
if ( $self->{_use_c}
&& $self->{_c_nodes}
&& $self->{c_psi} > 0
&& $threshold > 0
&& $threshold < 1 )
{
my $trees = $self->{trees};
my $t = scalar @$trees;
my $c = $self->{c_psi};
my ( $n_pts, $nf, $x_packed ) = $self->_resolve_input($data);
my $sums_packed = "\0" x ( $n_pts * 8 );
score_all_xs(
$self->{_c_nodes}, $self->{_c_coef_idx}, $self->{_c_coef_val},
$x_packed, $sums_packed, $n_pts,
$nf, $t, $self->{_use_openmp}
);
my $sum_threshold = -log($threshold) * $c * $t / log(2);
my $result = [];
predict_sums_xs( $sums_packed, $n_pts, $sum_threshold, $result );
return $result;
} ## end if ( $self->{_use_c} && $self->{_c_nodes} ...)
lib/Algorithm/Classifier/IsolationForest.pm view on Meta::CPAN
my ( $self, $data ) = @_;
$self->_check_fitted;
croak "pack_data requires the Inline::C backend; install Inline::C"
unless $self->{_use_c};
croak "pack_data() expects an arrayref of samples"
unless ref $data eq 'ARRAY';
my $n_pts = scalar @$data;
my $nf = $self->{n_features};
my $x_packed = "\0" x ( $n_pts * $nf * 8 );
my ( $mode, $fill ) = $self->_pack_args;
pack_input_xs( $data, $x_packed, $n_pts, $nf, $mode, $fill );
return bless {
packed => $x_packed,
n_pts => $n_pts,
n_feats => $nf,
},
'Algorithm::Classifier::IsolationForest::PackedData';
} ## end sub pack_data
# Internal helper: given $data that may be a raw arrayref OR a PackedData
# instance, return the (n_pts, n_feats, x_packed) triple ready for
# score_all_xs. Called from every scoring fast path.
sub _resolve_input {
my ( $self, $data ) = @_;
if ( ref $data eq 'Algorithm::Classifier::IsolationForest::PackedData' ) {
croak "PackedData has $data->{n_feats} features but model expects " . $self->{n_features}
unless $data->{n_feats} == $self->{n_features};
return ( $data->{n_pts}, $data->{n_feats}, $data->{packed} );
}
my $n_pts = scalar @$data;
my $nf = $self->{n_features};
my $x_packed = "\0" x ( $n_pts * $nf * 8 );
my ( $mode, $fill ) = $self->_pack_args;
pack_input_xs( $data, $x_packed, $n_pts, $nf, $mode, $fill );
return ( $n_pts, $nf, $x_packed );
} ## end sub _resolve_input
# Helper used by the pure-Perl fallback paths: convert either form back
# to an arrayref-of-arrayrefs. Slow on PackedData -- the whole point of
# packing is to keep things in C -- but lets the fallback path be
# uniformly arrayref-driven.
sub _to_arrayref {
my ( $self, $data ) = @_;
return $data if ref $data eq 'ARRAY';
if ( ref $data eq 'Algorithm::Classifier::IsolationForest::PackedData' ) {
my $n_pts = $data->{n_pts};
my $nf = $data->{n_feats};
my @doubles = unpack( 'd*', $data->{packed} );
my @rows;
for my $i ( 0 .. $n_pts - 1 ) {
push @rows, [ @doubles[ $i * $nf .. ( $i + 1 ) * $nf - 1 ] ];
}
return \@rows;
} ## end if ( ref $data eq 'Algorithm::Classifier::IsolationForest::PackedData')
croak "expected arrayref or PackedData, got " . ( ref($data) || 'scalar' );
} ## end sub _to_arrayref
# ---------------------------------------------------------------------------
# Missing-value handling.
#
# The `missing` strategy chosen at new() decides how undef feature cells are
# treated. Scoring always tolerates undef; the strategy governs fit() and
# how undef is represented for the scorer:
#
# die -- croak from fit() if the training data holds any undef cell.
# Scoring still maps undef -> 0 (the long-standing behaviour).
# zero -- undef counts as the value 0, at fit and score time.
# impute -- undef is replaced by a learned per-feature mean/median; the
# fill vector is stored on the model and reused at score time.
# nan -- ranges are built over present values only and a point missing
# the split feature is routed to the right child, consistently
# at fit (Perl) and score (C packs NaN; `<`/`<=` send it right).
# ---------------------------------------------------------------------------
# Returns the training data to actually build trees on, after applying the
# missing-value strategy. May croak (die), return a dense filled copy
# (zero/impute), or pass $data through unchanged (nan).
sub _prepare_fit_data {
my ( $self, $data ) = @_;
my $m = $self->{missing};
my $nf = $self->{n_features};
if ( $m eq 'die' ) {
for my $i ( 0 .. $#$data ) {
my $row = $data->[$i];
for my $f ( 0 .. $nf - 1 ) {
next if defined $row->[$f];
croak "fit(): undef feature value at sample $i, column $f; "
. "construct with missing => 'zero', 'impute', or 'nan' "
. "to train on data with missing values";
}
}
return $data;
} ## end if ( $m eq 'die' )
# nan: leave undef in place -- _build_tree / the split routers handle it.
return $data if $m eq 'nan';
# zero / impute: undef has to become a real number somewhere before a
# split can look at it. The fill vector is computed either way (it's
# needed for persistence and for scoring later), but densifying $data
# into a second, fully separate Perl array here is only necessary for
# the pure-Perl tree builder (_build_tree assumes every cell is
# defined once missing != 'nan' -- see its lo/hi scan). The C
# tree-building path -- _build_forest_c/_build_forest_openmp, and
# every parallel_fit worker, all of which go through pack_input_xs --
# already fills undef cells itself from this same fill vector, so
# skip the redundant whole-dataset copy when that's the path fit()
# will actually take. Scoring the training set for a learned
# contamination threshold (below, in fit()) is unaffected: it always
# runs through the pure-Perl scorer regardless of use_c (fit() drops
# any previous fit's packed buffers before that scoring, and
# _rebuild_c_trees runs after), and that path already tolerates raw
# undef cells
# for both zero (_path_length's "// 0") and impute (_prepare_perl_input
# densifies on demand from missing_fill).
my $fill
= $m eq 'impute'
? $self->_compute_impute_fill($data)
: [ (0) x $nf ];
$self->{missing_fill} = $fill if $m eq 'impute';
( run in 0.541 second using v1.01-cache-2.11-cpan-c966e8aa7e8 )