ベクトルと行列（linalg）

linalg パッケージは、ベクタと行列のための numpy スタイルの API を提供します。コンストラクタ、形状操作、要素ごとの算術演算、積、集約、離散微積分（差分と数値微分）、線形代数（行列式、逆行列、連立一次方程式の求解）が含まれます。

JSON ライブラリと同様に、linalg は Lisp ソース（linalg.lisp）として一度だけ実装されています。インタプリタは linalg: 関数が最初に使われたときに定義を遅延ロードし、コンパイルパスはパッケージが参照されたときに定義をプログラムへ継ぎ足します。バックエンドごとのコードは存在しないため、すべての関数はインタプリタ、JVM コンパイラ、WASM Preview 1、WASM コンポーネントで同一に振る舞います。

データ表現

linalg のコンストラクタは packed float 配列 を作ります。これは #d(...) リテラルと同じ、アンボックスな (array double-float) です。ベクタはランク 1 の配列で #d(1.0 2.0 ...) と印字され、行列はランク 2 の配列でネストした #d((...) ...) 形式で印字されます。#d はアンボックスな packed 表現を表すため、その印字結果を読み戻すと packed 配列になります。個々の要素は aref で読み書きでき、プログラムの他の場所で構築された配列 (packed でも一般のボックス配列でも) も linalg 関数に渡せます。より高いランクの配列も扱えます。要素ごとの演算、リダクション、reshape/flatten、array-equal はフラットな行優先順で要素を走査するため任意のランクを受け付けます。一方 dot/matmul/outer/det/inv/solve/trace/transpose は numpy の専用ルーチンと同様、ベクタと行列 (ランク 2 以下) に対して定義されたままです。linalg:from-list / linalg:to-list は配列とリストを相互に変換します。

linalg は速度を優先して浮動小数点で計算します。すべてのコンストラクタと配列を生成する演算はデフォルトで packed double-float 配列を返し（単精度も利用できます。単精度浮動小数点を参照）、linalg:det、linalg:inv、linalg:solve は (numpy と同様に) 浮動小数点で計算されるため、一般の逆行列には通常の丸めが生じ、ほぼ特異な行列式は厳密な 0 ではなく微小値になることがあります。リダクションは numpy と同様に要素の型に従います。packed または float 配列に対するリダクションは double を、素の整数配列 (#(1 2 3) のようなリテラル) に対するリダクションは整数または厳密な比を返します。linalg:norm は sqrt が浮動小数点数を返すため常に浮動小数点数です。クロスバックエンドの注意点が 1 つあります。WASM バックエンドは非終端の浮動小数点数をインタプリタや JVM より少ない有効桁数で印字するため、丸めのある逆行列や無理数のノルムは、内部の double が同一でもバックエンド間で見た目が異なることがあります。

実例

上記の inv と solve の行列は、float の結果が厳密になりどのバックエンドでも同一に印字されるよう選んでいます。(linalg:inv #2A((1 2) (3 4))) のような一般の逆行列は同じ値を計算しますが浮動小数点の丸めを伴います。

要素ごとの算術演算とブロードキャスト

linalg:add、linalg:sub、linalg:mul、linalg:div は要素ごとに演算し、numpy の規則でブロードキャストします: どちらか一方のスカラーのオペランドはもう一方のオペランドの形状にブロードキャストされ、形状の異なる 2 つの配列は末尾の軸から揃えられます -- 揃えた各軸の長さは等しいか、どちらかが 1 でなければならず（先頭側の欠けている軸は 1 と見なされます）、長さ 1 の軸がもう一方の長さに引き伸ばされます。どちらにも当てはまらないペアは shape mismatch エラーを通知します。結果は最初の配列オペランドの要素型を保持します（混合幅の規則と同じ）。mul はアダマール積（要素ごとの積）であることに注意してください。行列積は linalg:matmul（またはランクに応じてディスパッチする linalg:dot）です。要素ごとの任意の変換には linalg:emap を使います。

よく使う要素ごとの演算は、numpy の ufunc 名でも用意されています: linalg:exp・linalg:log・linalg:tanh・linalg:sin・linalg:cos・linalg:tan・linalg:asin・linalg:acos・linalg:atan・linalg:sinh・linalg:cosh・linalg:sqrt・linalg:abs・linalg:square・linalg:negative・linalg:sign・linalg:reciprocal、さらに比較セレクトの linalg:maximum・linalg:minimum・linalg:clip・linalg:relu(厳密比較 (if (> x y) x y) とその鏡像で定義され、比較が偽なら第 2 被演算子または境界が選ばれます — タイや NaN を含め、すべてのバックエンドで同一の規則です)。いずれも対応する emap(または mul / div / maximum / minimum の呼び出し)と等価ですが、名前付き関数なので --simd で加速されます。任意のコールバックを取る emap は決して加速されません。

軸に沿った還元

還元関数 linalg:sum、linalg:mean、linalg:amax、linalg:amin は省略可能な整数 axis(負は numpy 流に末尾から数える)を取り、配列全体ではなくその軸に沿って還元します。結果からその軸は除去されますが、省略可能な keepdims フラグが非 nil のときは長さ 1 の軸として保持されます — 入力にそのままブロードキャストで戻せる形状で、バッチ softmax が行ごとの最大値を引くのに使うのはこの形です。linalg:argmax と linalg:argmin も同じ axis 引数を取り、スライスごとの index を返します(行列に対しては packed double 配列 — linalg 配列に整数幅はありません)。linalg:reshape は 1 つの -1 extent を受け付け、要素数から推論します。

インデックス操作・選択・マスク

linalg:take-rows は index ベクタで axis-0 スライスを選択し(numpy の x[mask]、任意 rank)、linalg:gather は行ごとに 1 要素を取り出し(y[np.arange(n), t])、linalg:one-hot はラベル行列を作ります。要素ごとの比較 linalg:equal、linalg:greater、linalg:greater-equal、linalg:less、linalg:less-equal は 0.0/1.0 のマスクを返します(スカラー被演算子とブロードキャスト対応)。numpy が boolean index する場面では、マスクを掛け算してください。linalg:zeros-like は同じ形状・同じ幅のゼロ配列を確保します。

乱数

np.random に相当する乱数はシード可能で、バックエンド間で決定的です。linalg:seed は Wichmann-Hill 生成器をリセットし、その draw は正確な整数演算と IEEE double 演算だけで構成されるため、シード済みの linalg:rand、linalg:randn、linalg:uniform、linalg:choice、linalg:permutation の列は interpreter・JVM・両 WASM ターゲットで bit-identical です — 重み初期化とミニバッチ抽出がどこでも正確に再現されます。randn は Box-Muller ではなく Irwin-Hall(一様乱数 12 個の和)を使います(Box-Muller の log/cos は WASM で発散するため)。そのため裾は 6σ でクリップされます。初期化には十分ですが、分布が np.random.randn と厳密に一致するわけではありません。

離散微積分

linalg:diff と linalg:gradient は numpy の離散微積分ペア(np.diff / np.gradient)です。diff は最後の軸に沿った n 階の離散差分(デフォルト 1)を取ります。1 ステップごとにその軸が 1 つ短くなり、行列は各行内で差分されます。gradient はサンプル値のベクタの微分を、2 次精度の中心差分(両端は 1 次精度の片側差分)で推定するため、結果は入力と同じ長さになります。省略可能な第 2 引数には、一様なサンプル間隔(数値、デフォルト 1)か、非一様なサンプルのための同じ長さの座標ベクタを渡せます。どちらも他の linalg 変換と同様に入力の幅を保持します。算術は通常どおり浮動小数点ですが、厳密に微分できるサンプル値 — 以下の例のような、整数座標で読んだ多項式 — はすべてのバックエンドで同一に印字されます。

#(0 1 4 9 16) の gradient — 放物線 y = x^2 を x = 0..4 でサンプリングしたもの — は、内部の点で真の導関数 2x を厳密に復元します(中心差分は 2 次関数に対して厳密で、両端は 1 次精度の推定です)。座標ベクタ形式なら、最後の行のような不等間隔のサンプルでも厳密なままです。examples/ml/numerical-calculus.lisp は、これらの考え方を投射体の運動の例で一通り示します。

単精度浮動小数点

linalg はデフォルトで double-float で計算しますが、幅多相 (width-polymorphic) です。メモリが半分で SIMD レーン数が倍になる packed 単精度浮動小数点 (#f) 配列を受け付け、その幅を保持します。すべてのコンストラクタは末尾に省略可能な element-type を取り（デフォルトは 'double-float。#f の結果が欲しければ 'single-float を渡します）、すべての変換 -- add/sub/mul/div/emap、transpose/reshape、dot/matmul/outer、inv/solve -- は入力の幅を保持します。したがって単精度の値は最後まで単精度のまま流れます。関数的な重み更新 (linalg:sub W grad) は W の幅を保持し、暗黙に double へ戻す（JVM の --simd パスでは、続く vec:matvec で幅不一致エラーを強制する）ことはありません。f32 の速度とメモリが欲しく精度の低下を許容できるときは単精度を使い、det/inv/solve のような精度が重要な処理にはデフォルトの倍精度を使ってください。

SIMD アクセラレーション

linalg はどこでもフラグなしで正しく動きますが、--simd フラグで加速されます。32 の関数 — add・sub・mul・div・sum・norm・amax・amin・argmax・argmin・trace・transpose・reshape・dot・outer、単項 ufunc の exp・log・tanh・sin・cos・tan・asin・acos・atan・sinh・cosh・sqrt・abs・negative・sign、および比較セレクトの maximum・minimum — がネイティブなベクトルカーネル(インタプリタと JVM では jdk.incubator.vector、wasm-GC では WebAssembly の v128)にルーティングされ、それらを使って書かれている mean・matmul・flatten・solve・square・reciprocal・clip・relu も一緒に加速されます。プログラムが何を受け付け何を拒否するかは一切変わりません。カーネルが扱えない入力(一般の boxed 配列、幅の混在、素の数値)は移植可能な linalg.lisp の定義がそのまま実行され、同じ結果と同じエラーメッセージになります。観測可能な違いは単精度リダクションの精度規則だけで、要素ごとの結果と完全な行列積はビット一致のままです。

速度のためにパッケージを乗り換える理由はありません。--simd の下では linalg と vec は同じカーネルに行き着きます。vec と linalg の使い分けを参照してください — 要約すると: デフォルトでは linalg に対して書き、vec に手を伸ばすのは -into の書き込み先渡しループ・--no-gc ターゲット(linalg はそこではコンパイルできません)・fail-fast な幅の厳格さが必要なときだけです。

第一級関数

linalg の関数は通常の defun なので、#'linalg:norm などは関数が期待されるあらゆる場所で第一級の値として動作します。

配列は同一性 (eq) でしか比較されないため、結果の比較には形状と数値の等価性を検査する linalg:array-equal を使います（1 と 1.0 は等しいと判定されます）。

パッケージ

linalg はそれ自体が独立したパッケージであり、cl を use していません。(in-package linalg) の中では標準関数は非修飾の名前では見えず、cl: 修飾（cl:print、cl:mapcar など）が必要になります。したがってほとんどのプログラムは、このページのすべての例がそうしているように、デフォルトの cl-user パッケージにとどまり、修飾された linalg: 名で呼び出すべきです。