Arraymancer icon indicating copy to clipboard operation
Arraymancer copied to clipboard

Export `cumsum` procedure and expand to work on n-dim Tensors

Open filipeclduarte opened this issue 3 years ago • 1 comments

I have three suggestions:

  1. Export cumsum procedure: https://github.com/mratsim/Arraymancer/blob/1a2422a1e150a9794bfaa28c62ed73e3c7c41e47/src/arraymancer/ml/clustering/kmeans.nim#L26
proc cumsum[T: SomeFloat](p: Tensor[T]): Tensor[T] {.noInit.} =
  ## Calculates the cumulative sum of a vector.
  ## Inputs:
  ##  - p: a rank-1 tensor to cumulatively sum
  ##
  ## TODO: implement parallel prefix sum
  ## See: https://en.wikipedia.org/wiki/Prefix_sum#Algorithm_2:_Work-efficient
  ##
  ## Returns:
  ##  - A tensor cumulatively summed, that is, add each value to
  ##    all previous values, sequentially
  result = p.clone()
  assert p.rank == 1
  assert p.shape[1] == 0
  let n_rows = p.shape[0]
  for i in 1..<n_rows:
    result[i] += result[i-1]
  1. Generalize this procedure to work on n-dimensional Tensors.

  2. Create a procedure to do cumulative product cumprod that is similar to cumsum, but instead of sum, use prod.

filipeclduarte avatar Jan 04 '21 23:01 filipeclduarte

@HugoGranstrom suggested this proc cumprod at the discord:

proc cumsum*[T](t: Tensor[T], axis: int): Tensor[T] =
  ## Calculates the cumulative sum of a rank-n Tensor.
  ## Inputs:
  ##  - t: a rank-n tensor to cumulatively sum
  ##  - axis: int
  ## Returns:
  ##  - A tensor cumulatively summed at axis
  result = zeros_like(t)
  for i, tAxis in enumerateAxis(t, axis):
    var temp = result.atAxisIndex(axis, i)
    if i == 0:
      temp[_] = tAxis
    else:
      temp[_] = result.atAxisIndex(axis, i-1) + tAxis

Then I suggested to create the cumprod as:

proc cumprod*[T](t: Tensor[T], axis:int): Tensor[T] = # from hugogranstrom
  ## Calculates the cumulative product of a rank-n Tensor.
  ## Inputs:
  ##  - t: a rank-n tensor to cumulatively sum
  ##  - axis: int
  ## Returns:
  ##  - A tensor cumulatively product at axis, that is, product each value
  result = zeros_like(t)
  for i, tAxis in enumerateAxis(t, axis):
    var temp = result.atAxisIndex(axis, i)
    if i == 0:
      temp[_] = tAxis
    else:
      temp[_] = result.atAxisIndex(axis, i-1) *. tAxis

I think these are useful for numerous tasks in scientific computing and finance.

filipeclduarte avatar Jan 05 '21 01:01 filipeclduarte