Class: CArray::Histogram

Inherits:
Object
  • Object
show all
Defined in:
lib/carray/histogram.rb

Overview

Binned counts over one or more continuous variables, of any dimensionality M. Built by CArray#histogram1d / #histogram2d / #histogram rather than constructed directly; the 1-D entry point is the same class with M = 1.

For discrete integer labels (value == bin index, no edges) use BincountND instead.

Instance Attribute Summary collapse

Instance Method Summary collapse

Constructor Details

#initialize(edges:, fiber_shape: [], include_max: false, weights_dtype: nil) ⇒ Histogram

Returns a new instance of Histogram.

Allocates a new histogram accumulator.

Parameters:

  • edges (Array<CArray, Array<Numeric>>) —

    one edges array per histogram dimension; each must be 1-D sorted ascending with at least 2 values.

  • fiber_shape (Array<Integer>) (defaults to: []) —

    shape of the leading (non-sample, non-channel) axes.

  • include_max (Boolean, Array<Boolean>) (defaults to: false) —

    whether values equal to the last edge fold into the last bin; a scalar broadcasts across dimensions.

  • weights_dtype (Symbol, nil) (defaults to: nil) —

    data_type of the accumulator when weighted; nil for a count-only accumulator (int64 counts).

Raises:

  • (ArgumentError) 


152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
 
# File 'lib/carray/histogram.rb', line 152

def initialize (edges:, fiber_shape: [], include_max: false, weights_dtype: nil)
  @edges_list = edges.map { |e| CArray.wrap_readonly(e, :float64) }
  raise ArgumentError, "edges must be a non-empty list" if @edges_list.empty?
  @edges_list.each_with_index do |e, k|
    raise ArgumentError, "edges[#{k}] must be 1-D" unless e.ndim == 1
    raise ArgumentError, "edges[#{k}] needs at least 2 values" if e.elements < 2
  end
  @m = @edges_list.size                     # histogram dimensionality (= channel axis length)
  @n_list = @edges_list.map { |e| e.elements - 1 }   # per-dim bin count
  @fiber_shape = fiber_shape.map(&:to_i).freeze
  @include_max = case include_max
                 when Array
                   raise ArgumentError, "include_max length mismatch" unless include_max.size == @m
                   include_max.map { |v| !!v }
                 else
                   [!!include_max] * @m
                 end
  @weighted = !weights_dtype.nil?
  @counts_dtype = @weighted ? weights_dtype : :int64
  ext_dims = @n_list.map { |n| n + 2 }
  ext_shape = @fiber_shape + ext_dims
  @full_counts = CArray.public_send(@counts_dtype, *ext_shape).fill(0)
  @sample_axis  = nil
  @channel_axis = nil
end

Instance Attribute Details

#edges_list ⇒ Object (readonly)

Returns the value of attribute edges_list.



179
180
181
 
# File 'lib/carray/histogram.rb', line 179

def edges_list
  @edges_list
end

#fiber_shape ⇒ Object (readonly)

Returns the value of attribute fiber_shape.



179
180
181
 
# File 'lib/carray/histogram.rb', line 179

def fiber_shape
  @fiber_shape
end

#full_counts ⇒ Object (readonly)

Returns the value of attribute full_counts.



179
180
181
 
# File 'lib/carray/histogram.rb', line 179

def full_counts
  @full_counts
end

#include_max ⇒ Object (readonly)

Returns the value of attribute include_max.



179
180
181
 
# File 'lib/carray/histogram.rb', line 179

def include_max
  @include_max
end

#m ⇒ Object (readonly)

Returns the value of attribute m.



179
180
181
 
# File 'lib/carray/histogram.rb', line 179

def m
  @m
end

#n_list ⇒ Object (readonly)

Returns the value of attribute n_list.



179
180
181
 
# File 'lib/carray/histogram.rb', line 179

def n_list
  @n_list
end

Instance Method Details

#+(other) ⇒ Histogram

Returns a new Histogram whose counts are the element-wise sum of self and other. Both operands must share edges, fiber_shape, include_max, and weighted/unweighted state.

Parameters:

  • other (Histogram) —

    compatible accumulator.

Returns:

Raises:

  • (ArgumentError) —

    when the structure does not match.

Raises:

  • (ArgumentError) 


367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
 
# File 'lib/carray/histogram.rb', line 367

def + (other)
  # --- semantic guards: structure must match exactly ---------------
  raise ArgumentError, "type mismatch" unless other.is_a?(Histogram)
  raise ArgumentError, "M mismatch" unless @m == other.m
  @edges_list.each_with_index do |e, k|
    raise ArgumentError, "edges[#{k}] mismatch" unless e == other.edges_list[k]
  end
  raise ArgumentError, "fiber_shape mismatch" unless @fiber_shape == other.fiber_shape
  raise ArgumentError, "include_max mismatch (semantic guard)" unless @include_max == other.include_max
  raise ArgumentError, "weighted/unweighted mismatch" unless @weighted == other.weighted?

  result = self.class.send(:new,
                           edges: @edges_list,
                           fiber_shape: @fiber_shape,
                           include_max: @include_max,
                           weights_dtype: @weighted ? @counts_dtype : nil)
  rf = result.instance_variable_get(:@full_counts)
  rf[] = @full_counts + other.full_counts
  result.instance_variable_set(:@sample_axis, @sample_axis)
  result.instance_variable_set(:@channel_axis, @channel_axis)
  result
end

#add(chunk, axis: nil, weights: nil) ⇒ self

Accumulates chunk into self. On the first call the sample and channel axes are locked; subsequent calls must supply the same axis pair. When the accumulator is weighted, weights are required with a shape equal to chunk.shape minus the channel axis.

Parameters:

  • chunk (CArray) —

    sample values with shape fiber_shape + (A, M) (channel axis size must equal m).

  • axis (Array(Integer, Integer), Integer, nil) (defaults to: nil) —

    [sample, channel] axis pair; a bare Integer is treated as the sample axis for 1-D accumulators.

  • weights (CArray, nil) (defaults to: nil) —

    optional per-sample weights.

Returns:

  • (self)

Raises:

  • (ArgumentError) —

    on shape / axis / weighted-state mismatch.

Raises:

  • (ArgumentError) 


266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
 
# File 'lib/carray/histogram.rb', line 266

def add (chunk, axis: nil, weights: nil)
  chunk = CArray.wrap_readonly(chunk, :float64)

  # For M=1, accept chunks without the trailing channel axis (= the 1-D
  # user convention from `histogram1d`).  Auto-reshape adds a length-1
  # axis at the end; scalar `axis:` is interpreted as the sample axis in
  # the unwrapped layout.
  if @m == 1 && chunk.ndim == @fiber_shape.size + 1
    chunk = chunk.reshape(*(chunk.shape + [1]))
    if axis.is_a?(Integer)
      # axis was given in the unwrapped (pre-reshape) layout: normalize
      # against ndim-1 (= the unwrapped ndim) then pair with the new
      # trailing channel position.
      ax = CArray.normalize_axis(axis, chunk.ndim - 1, "add axis")
      axis = [ax, chunk.ndim - 1]
    end
  end

  # --- normalize axis: into [sample, channel] pair -----------------
  ax = axis || [-2, -1]
  ax = [ax] if ax.is_a?(Integer)
  raise ArgumentError, "axis must be [sample, channel]" unless ax.is_a?(Array) && ax.size == 2
  sample_ax  = CArray.normalize_axis(ax[0], chunk.ndim, "sample axis")
  channel_ax = CArray.normalize_axis(ax[1], chunk.ndim, "channel axis")
  raise ArgumentError, "same axis used twice" if sample_ax == channel_ax

  # --- lock axes on first add, otherwise verify against locked -----
  if @sample_axis.nil?
    @sample_axis = sample_ax
    @channel_axis = channel_ax
  elsif @sample_axis != sample_ax || @channel_axis != channel_ax
    raise ArgumentError,
          "axis mismatch (locked at [#{@sample_axis}, #{@channel_axis}], got [#{sample_ax}, #{channel_ax}])"
  end

  # --- validate chunk shape against (fiber_shape, M) ---------------
  expected_ndim = @fiber_shape.size + 2
  unless chunk.ndim == expected_ndim
    raise ArgumentError,
          "chunk.ndim=#{chunk.ndim} expected #{expected_ndim} " \
          "(fiber #{@fiber_shape.inspect} + sample + channel)"
  end
  unless chunk.shape[channel_ax] == @m
    raise ArgumentError,
          "channel axis length #{chunk.shape[channel_ax]} != M=#{@m}"
  end
  chunk_fiber = chunk.shape.dup
  [sample_ax, channel_ax].sort.reverse.each { |p| chunk_fiber.delete_at(p) }
  unless chunk_fiber == @fiber_shape
    raise ArgumentError,
          "fiber shape mismatch: chunk yields #{chunk_fiber.inspect}, expected #{@fiber_shape.inspect}"
  end

  sample_count = chunk.shape[sample_ax]
  return self if sample_count == 0

  if weights
    raise ArgumentError, "weights given but accumulator is unweighted" unless @weighted
    weights = CArray.wrap_readonly(weights, @counts_dtype)
    expected_w_shape = chunk.shape.dup
    expected_w_shape.delete_at(channel_ax)
    unless weights.shape == expected_w_shape
      raise ArgumentError,
            "weights shape #{weights.shape.inspect} expected #{expected_w_shape.inspect} " \
            "(chunk minus channel axis at #{channel_ax})"
    end
  elsif @weighted
    raise ArgumentError, "weights required (accumulator is weighted)"
  end

  # --- fused scatter kernel (stage 2) ------------------------------
  # Bin all M channels per sample and scatter directly into @full_counts
  # with NO intermediate index arrays (peak memory O(1), not O(M * A) for
  # A samples).
  # self is transposed to [fiber..., sample, channel] (a view; the kernel
  # iterator delivers it strided, no materialise).  Weights, if present,
  # are transposed to [fiber..., sample] and delivered by a second
  # iterator in lockstep (both views, no materialise).
  fiber_axes = (0...chunk.ndim).to_a - [sample_ax, channel_ax]
  tchunk = chunk.transpose(*(fiber_axes + [sample_ax, channel_ax]))

  tweights = nil
  if weights
    # weights axes = chunk axes with channel removed: an index above
    # channel_ax shifts down by 1.
    shift = ->(p) { p < channel_ax ? p : p - 1 }
    tweights = weights.transpose(*(fiber_axes.map(&shift) + [shift.call(sample_ax)]))
  end

  tchunk.send(:histogram_scatter_ki, @full_counts, @edges_list, @include_max, tweights)

  self
end

#counts ⇒ CArray

Returns the in-range counts view with shape fiber_shape + (n1, n2, ..., nM), excluding under- and over-flow bins.

Returns:



194
195
196
197
 
# File 'lib/carray/histogram.rb', line 194

def counts
  idx = [nil] * @fiber_shape.size + [1..-2] * @m
  @full_counts[*idx]
end

#edges ⇒ CArray+

Returns the bin edges: a single CArray when the accumulator is 1-D (M == 1), an Array of CArrays otherwise.

Returns:



185
186
187
 
# File 'lib/carray/histogram.rb', line 185

def edges
  @m == 1 ? @edges_list[0] : @edges_list
end

#midpoints ⇒ CArray+

Returns the midpoint of each in-range bin. Polymorphic like #edges: a single CArray for M == 1, an Array of CArrays otherwise.

Returns:



229
230
231
232
 
# File 'lib/carray/histogram.rb', line 229

def midpoints
  arr = @edges_list.map { |e| (e[0..-2] + e[1..-1]) / 2.0 }
  @m == 1 ? arr[0] : arr
end

#outlier_total ⇒ CArray

Returns the per-fiber count of samples that fell outside every in-range bin.

Returns:



247
248
249
 
# File 'lib/carray/histogram.rb', line 247

def outlier_total
  sum_along_bin_axes(@full_counts) - sum_along_bin_axes(counts)
end

#over(axis: nil) ⇒ CArray

Returns the overflow marginal along the given bin axis, shape fiber_shape. Same axis convention as #under.

Parameters:

  • axis (Integer, nil) (defaults to: nil) —

    bin dimension to marginalise.

Returns:

Raises:

  • (ArgumentError) —

    when axis is required but omitted.

Raises:

  • (ArgumentError) 


218
219
220
221
222
 
# File 'lib/carray/histogram.rb', line 218

def over (axis: nil)
  raise ArgumentError, "axis: keyword required (M=#{@m})" if axis.nil? && @m > 1
  ax = axis.nil? ? 0 : CArray.normalize_axis(axis, @m, "over")
  outlier_marginal(ax, -1)
end

#total ⇒ CArray

Returns the per-fiber sample total, including outliers, with shape fiber_shape (or a scalar when fiber_shape is empty).

Returns:



239
240
241
 
# File 'lib/carray/histogram.rb', line 239

def total
  sum_along_bin_axes(@full_counts)
end

#under(axis: nil) ⇒ CArray

Returns the underflow marginal along the given bin axis, shape fiber_shape. For 1-D accumulators axis may be omitted; for joint histograms it must be specified.

Parameters:

  • axis (Integer, nil) (defaults to: nil) —

    bin dimension to marginalise.

Returns:

Raises:

  • (ArgumentError) —

    when axis is required but omitted.

Raises:

  • (ArgumentError) 


206
207
208
209
210
 
# File 'lib/carray/histogram.rb', line 206

def under (axis: nil)
  raise ArgumentError, "axis: keyword required (M=#{@m})" if axis.nil? && @m > 1
  ax = axis.nil? ? 0 : CArray.normalize_axis(axis, @m, "under")
  outlier_marginal(ax, 0)
end