エンジニアのソフトウェア的愛情

または私は如何にして心配するのを止めてプログラムを・愛する・ようになったか

図形で理解するということ〜線形代数と数的霊薬〜

行列の積

 3 \times 3 の行列同士の積を考えます。


\begin{bmatrix}
a_{11} & a_{12} & a_{13} \\
a_{21} & a_{22} & a_{23} \\
a_{31} & a_{32} & a_{33}
\end{bmatrix}

\cdot

\begin{bmatrix}
b_{11} & b_{12} & b_{13} \\
b_{21} & b_{22} & b_{23} \\
b_{31} & b_{32} & b_{33}
\end{bmatrix}

定義に従い計算すると、1行1列の要素は左の項の1行目と右の項の1列目の各々の要素の積の和、1行2列の要素は左の項の1行目と右の項の2列目の各々の要素の積の和、…、3行3列の要素は左の項の3行目と右の項の3列目の各々の要素の積の和、となります。


\begin{bmatrix}
a_{11} b_{11} + a_{12} b_{21} + a_{13} b_{31} & a_{11} b_{12} + a_{12} b_{22} + a_{13} b_{32} & a_{11} b_{13} + a_{12} b_{23} + a_{13} b_{33} \\
a_{21} b_{11} + a_{22} b_{21} + a_{23} b_{31} & a_{21} b_{12} + a_{22} b_{22} + a_{23} b_{32} & a_{21} b_{13} + a_{22} b_{23} + a_{23} b_{33} \\
a_{31} b_{11} + a_{32} b_{21} + a_{33} b_{31} & a_{31} b_{12} + a_{32} b_{22} + a_{33} b_{32} & a_{31} b_{13} + a_{32} b_{23} + a_{33} b_{33}
\end{bmatrix}

この行列の和の部分を行列の外に出すと、次のように表現することができます。


\begin{bmatrix}
a_{11} b_{11} & a_{11} b_{12} & a_{11} b_{13} \\
a_{21} b_{11} & a_{21} b_{12} & a_{21} b_{13} \\
a_{31} b_{11} & a_{31} b_{12} & a_{31} b_{13} 
\end{bmatrix}

+

\begin{bmatrix}
a_{12} b_{21} & a_{12} b_{22} & a_{12} b_{23} \\
a_{22} b_{21} & a_{22} b_{22} & a_{22} b_{23} \\
a_{32} b_{21} & a_{32} b_{22} & a_{32} b_{23}
\end{bmatrix}

+

\begin{bmatrix}
a_{13} b_{31} & a_{13} b_{32} & a_{13} b_{33} \\
a_{23} b_{31} & a_{23} b_{32} & a_{23} b_{33} \\
a_{33} b_{31} & a_{33} b_{32} & a_{33} b_{33}
\end{bmatrix}

これは、左の行列を列ベクトルに、右の行列を行ベクトルに分けたとき、1 番目の列ベクトルと行ベクトル、 2 番目の列ベクトルと行ベクトル、3 番目の列ベクトルと行ベクトルのそれぞれの積を足し合わせたものと同じになります。


\begin{array}{rcl}

\begin{bmatrix}
| & | & | \\
\boldsymbol{a_{1}} & \boldsymbol{a_{2}} & \boldsymbol{a_{3}} \\
| & | & |
\end{bmatrix}

\cdot

\begin{bmatrix}
- & \boldsymbol{b^*_{1}} & - \\
- & \boldsymbol{b^*_{2}} & - \\
- & \boldsymbol{b^*_{3}} & -
\end{bmatrix}

& = &
\boldsymbol{a_1} \boldsymbol{b^*_1} +
\boldsymbol{a_2} \boldsymbol{b^*_2} +
\boldsymbol{a_3} \boldsymbol{b^*_3}

\end{array}

これを計算で確かめます。

Nx - Numerical Elixir

Elixir には Numerical Elixir - Nx というプロジェクト / ライブラリがあります。

github.com

nx.hexdocs.pm

機械学習などで利用されることを想定されていますが、ここでは純粋に行列計算に利用します。

まず IEx を起動します。

$ iex
Erlang/OTP 29 [erts-17.0.5] [source] [64-bit] [smp:8:8] [ds:8:8:10] [async-threads:1] [jit] [dtrace]

Interactive Elixir (1.20.4) - press Ctrl+C to exit (type h() ENTER for help)
iex(1)> 

IEx が起動したら、環境に Nx をインストールします。

iex> Mix.install([:nx])

以下、プロンプトの iex> は省略します。

ここでは Nx モジュールの関数をモジュール名の装飾なしに利用したいので、すべて import します。

import Nx

行列を定義します。

a = ~MAT[
1 2 3
4 5 6
7 8 9
]
#Nx.Tensor<
  s32[3][3]
  [
    [1, 2, 3],
    [4, 5, 6],
    [7, 8, 9]
  ]
>

ここで ~MATNx.Tensor の値を作成するための構文糖衣です。 あとで出てくる ~VEC も同様です。

行列をもう一つ定義します。

b = ~MAT[
1 4 7
2 5 8
3 6 9
]
#Nx.Tensor<
  s32[3][3]
  [
    [1, 4, 7],
    [2, 5, 8],
    [3, 6, 9]
  ]
>

行列から列ベクトルを取り出します。

take(a, ~VEC[0], axis: 1)
#Nx.Tensor<
  s32[3][1]
  [
    [1],
    [4],
    [7]
  ]
>

axis には軸を指定します。 表示された Nx.Tensor の値からわかる通り、値は [[1,2,3],[4,5,6],[7,8,9]] という並びで格納されています。 Nx では添え字は 0 から始まるので、縦方向に値を拾いたい場合は axis: 1 を指定します。

同様に行ベクトルを取り出します。 axis オプションのデフォルト値は 0 のため指定を省略できますが、ここでは表記を合わせることにします。

take(b, ~VEC[0], axis: 0)
#Nx.Tensor<
  s32[1][3]
  [
    [1, 4, 7]
  ]
>

なお第 2 引数にはスカラも利用できますが、ここでスカラを与えると期待する形の値が得られないため、それぞれベクトル( ~VEC[0] )で指定しています。

take(b, 0, axis: 0)
#Nx.Tensor<
  s32[3]
  [1, 4, 7]
>

行列の積の値を確認します。

dot(a, b)
#Nx.Tensor<
  s32[3][3]
  [
    [14, 32, 50],
    [32, 77, 122],
    [50, 122, 194]
  ]
>

次に、各々の列ベクトルと行ベクトルの積の和の値を確認します。

a1b1 = dot(take(a, ~VEC[0], axis: 1), take(b, ~VEC[0], axis: 0))
a2b2 = dot(take(a, ~VEC[1], axis: 1), take(b, ~VEC[1], axis: 0))
a3b3 = dot(take(a, ~VEC[2], axis: 1), take(b, ~VEC[2], axis: 0))

a1b1 |> add(a2b2) |> add(a3b3)
#Nx.Tensor<
  s32[3][3]
  [
    [14, 32, 50],
    [32, 77, 122],
    [50, 122, 194]
  ]
>

一致することが確認できました。

Julia - The Julia Programming Language

同じ内容を Julia でも計算させてみます。

julialang.org

$ julia
               _
   _       _ _(_)_     |  Documentation: https://docs.julialang.org
  (_)     | (_) (_)    |
   _ _   _| |_  __ _   |  Type "?" for help, "]?" for Pkg help.
  | | | | | | |/ _` |  |
  | | |_| | | | (_| |  |  Version 1.7.3 (2022-05-06)
 _/ |\__'_|_|_|\__'_|  |  Official https://julialang.org/ release
|__/                   |

以下、プロンプトの julia> 以降が入力した式で、その下に表示される値は Julia の出力です。

julia> A = [1 2 3; 4 5 6; 7 8 9]
3×3 Matrix{Int64}:
 1  2  3
 4  5  6
 7  8  9

julia> B = [1 4 7; 2 5 8; 3 6 9]
3×3 Matrix{Int64}:
 1  4  7
 2  5  8
 3  6  9
julia> A * B
3×3 Matrix{Int64}:
 14   32   50
 32   77  122
 50  122  194

Julia では、ベクトルの後に [行番号,列番号] と列番号と列番号を指定することで、指定した位置の値を取得できます。 値の代わりに : を指定すると列全体または行全体を取得できます。

なお Julia では添え字は 1 始まりなので、1 列目や 1 行目には 1 を指定します。

ただし、単純に列番号を指定すると扱いたい形式にならないため、

julia> A[1,:]
3-element Vector{Int64}:
 1
 2
 3

[[1],:] という形で指定します。

julia> A[[1],:]
1×3 Matrix{Int64}:
 1  2  3
julia> A[:,1]
3-element Vector{Int64}:
 1
 4
 7

取り出した列と行は、同じように積をとることができます。

julia> A[:,1] * B[[1],:]
3×3 Matrix{Int64}:
 1   4   7
 4  16  28
 7  28  49

各々の列ベクトルと行ベクトルの積の和の値を確認します。

julia> A[:,1] * B[[1],:] + A[:,2] * B[[2],:] + A[:,3] * B[[3],:]
3×3 Matrix{Int64}:
 14   32   50
 32   77  122
 50  122  194

一致することが確認できました。

図解線形代数

なお、ちゃんとした話はこちらにありますので、こちらを読んでいただけたらと。

anagileway.com