Luna Flow
日本語

algebra チュートリアル

このチュートリアルでは、線形代数のヘルパーを algebra の trait に対して一度だけ書き、その trait を実装する任意の行列型・ベクトル型（リポジトリの密なラッパーや自作の型）で実行する方法を示します。MoonBit の総称型の知識が必要です。数学は軽めにとどめ、algebra の設計 で説明しています。

クイックスタート

モジュールと、使用するパッケージを追加します。

moon add Luna-Flow/linear-algebra@0.5.0
///|
import {
  "Luna-Flow/linear-algebra/algebra",
  "Luna-Flow/linear-algebra/backends/default",
}

最小限の有用なプログラムは、trait を通じて行列に形状を問い合わせます。そのため、どの MatrixShape 型でも動作します。

///|
fn[M : @algebra.MatrixShape] alg_tut_describe(matrix : M) -> String {
  let (rows, cols) = @algebra.MatrixShape::shape(matrix)
  "\{rows}x\{cols}"
}

///|
test "describe a dense matrix" {
  let m = @default.DenseMatrix::from_2d_array([[1, 2, 3], [4, 5, 6]])
  inspect(alg_tut_describe(m), content="2x3")
}

出力は 2x3 です。

日常的なタスク

任意の加法的ベクトル型で残差を計算する

残差 r=b−Axr = b - A x に必要なのはベクトルの減算だけです。AdditiveVector を要求すれば、ヘルパーは可変と不変の両方の密ベクトルを受け付けます。

///|
fn[V : @algebra.AdditiveVector] alg_tut_residual(
  observed : V,
  predicted : V,
) -> V {
  observed - predicted
}

///|
test "residual over two vector representations" {
  let b = @default.DenseVector::from_array([3, 5, 7])
  let p = @default.DenseVector::from_array([1, 5, 9])
  let r = alg_tut_residual(b, p)
  inspect(r.inner(), content="|2, 0, -2|")
  let bi = @default.ImmutableDenseVector::from_array([3, 5, 7])
  let pi = @default.ImmutableDenseVector::from_array([1, 5, 9])
  inspect(alg_tut_residual(bi, pi).inner(), content="|2, 0, -2|")
}

Gram 行列 ATAA^{\mathsf T} A を構築する

Gram 行列には転置と行列積が必要で、それがまさに MatMulMatrix です。

///|
fn[M : @algebra.MatMulMatrix] alg_tut_gram(a : M) -> M {
  @algebra.TransposeMatrix::transpose(a) * a
}

///|
test "gram matrix of a 3x2 design matrix" {
  let a = @default.ImmutableDenseMatrix::from_2d_array([[1, 0], [1, 1], [1, 2]])
  let g = alg_tut_gram(a)
  debug_inspect(@algebra.MatrixShape::shape(g), content="(2, 2)")
  inspect(g.inner(), content="|3, 3|\n|3, 5|")
}

ATAA^{\mathsf T}A は AA の形状にかかわらず常に定義されます。m×nm \times n 行列からは n×nn \times n の結果が得られます。

乗算の前に合成可能性を検査する

実行時に形状が決まる行列では * は部分的であり、密なラッパーは不一致で中断します。形状がデータから来る場合は、先に MatrixShape で検査してください。

///|
fn[M : @algebra.MatMulMatrix] alg_tut_try_product(a : M, b : M) -> M? {
  let (_, inner_a) = @algebra.MatrixShape::shape(a)
  let (inner_b, _) = @algebra.MatrixShape::shape(b)
  if inner_a == inner_b {
    Some(a * b)
  } else {
    None
  }
}

///|
test "product only for composable shapes" {
  let a = @default.DenseMatrix::from_2d_array([[1, 2]])
  let b = @default.DenseMatrix::from_2d_array([[3], [4]])
  debug_inspect(
    alg_tut_try_product(a, b).map(m => m.inner().get(0, 0)),
    content="Some(11)",
  )
  inspect(alg_tut_try_product(a, a) is None, content="true")
}

自作の固定サイズ行列を使う

2×22 \times 2 行列型は積が全域的なので、実行時の前提条件なしにすべての行列の層を実装できます。演算子 trait を実装してから、各層を宣言します。

///|
struct AlgTutMat2 {
  a : Int
  b : Int
  c : Int
  d : Int
}

///|
impl @algebra.MatrixShape for AlgTutMat2 with fn shape(_) {
  (2, 2)
}

///|
impl @algebra.TransposeMatrix for AlgTutMat2 with fn transpose(m) {
  { a: m.a, b: m.c, c: m.b, d: m.d, }
}

///|
impl Add for AlgTutMat2 with fn add(x, y) {
  { a: x.a + y.a, b: x.b + y.b, c: x.c + y.c, d: x.d + y.d, }
}

///|
impl Neg for AlgTutMat2 with fn neg(x) {
  { a: -x.a, b: -x.b, c: -x.c, d: -x.d, }
}

///|
impl Sub for AlgTutMat2 with fn sub(x, y) {
  x + -y
}

///|
impl Mul for AlgTutMat2 with fn mul(x, y) {
  {
    a: x.a * y.a + x.b * y.c,
    b: x.a * y.b + x.b * y.d,
    c: x.c * y.a + x.d * y.c,
    d: x.c * y.b + x.d * y.d,
  }
}

///|
impl @algebra.AdditiveMatrix for AlgTutMat2

///|
impl @algebra.MatMulMatrix for AlgTutMat2

///|
test "the generic gram helper runs on a custom type" {
  let m : AlgTutMat2 = { a: 1, b: 2, c: 3, d: 4, }
  let g = alg_tut_gram(m)
  debug_inspect((g.a, g.b, g.c, g.d), content="(10, 14, 14, 20)")
}

前のタスクのヘルパー alg_tut_gram は変更不要でした。

さらに進んで

backends/default の汎用ヘルパーを使う。 @default.shape_of、@default.transpose、@default.matmul は、3 つの演算の trait 境界付き版です。trait 修飾の呼び出しではなく関数値が欲しいときに便利です。

container と組み合わせる。 algebra の trait はオブジェクト全体の操作を記述し、要素を公開することはありません。アルゴリズムが個々の要素の読み取りや構築も必要とする場合は、trait 境界と並べて container の操作辞書を使ってください。2 つの層は独立しています。

スカラーの要件は具体的な型にとどまる。 trait はスカラー型に言及しません。アルゴリズムがスカラー倍や内積を必要とする場合は、具体的な型（T : AddMonoid + Mul を伴う @default.DenseVector[T]）を受け取るか、スカラー演算を関数として渡してください。

検査付きの積。 失敗を値にしなければならない場合は、Result[_, LinearAlgebraError] を返す @immut.Matrix::matmul のような検査付きメソッドを持つ具体的な型に変換してください。error のチュートリアル を参照してください。

よくある落とし穴

  • 汎用コードで trait メソッドをドット構文で呼ぶ。 fn[M : @algebra.MatMulMatrix] の中では @algebra.TransposeMatrix::transpose(m) と書いてください。型パラメータに対してスーパートレイトのメソッドをドット呼び出しすることは、MoonBit 0.10 で非推奨になっています。
  • Hadamard 積の * で MatMulMatrix を実装する。 MatMulMatrix の * は行列積でなければなりません。* が要素ごとの型は AdditiveMatrix までにとどめてください。
  • あらゆるスカラーで (AB)^T = B^T A^T を期待する。 これには可換なスカラーが必要です。四元数行列ではこれが成り立ちません。
  • 浮動小数点の結果を厳密に比較する。 Double 行列の積は丸め誤差の範囲でしか一致しません。許容誤差付きで比較してください（arithmetic のチュートリアル を参照）。

次のステップ