moGue’s Garten

わりとアレな知識をアレな人向けになりふり構わず放出する

numpy の気づき

numpy では,こういうこともできる.

import numpy as np

A = np.array([[1, 2, -3, 4], 
              [5, -6, 7, 8], 
              [9, 10, 11, -12]])

A_sq = A ** 2
A_abs = np.abs(A)

A_positive = (A>0)
#[[ True  True False  True]
# [ True False  True  True]
# [ True  True  True False]]

X = np.empty_like(A)
X[A_positive] = A_sq[A_positive]
X[~A_positive] = A_abs[~A_positive]

print(X)
# [[  1   4   3  16]
#  [ 25   6  49  64]
#  [ 81 100 121  12]]

こういうのはスライスというのか. 配列の要素ごとに演算を行いたいが, 演算に場合分けを要する場合などに使えそうだ. ただし,A_sq および A_abs によって結局2通りの場合を全要素に対して試すことになるから,計算の無駄はある. それでも,たいていの計算ならば for 文を使うよりは速そう.

matplotlib を学ぶ

python のグラフ描画ツールといえば matplotlib. python で計算してえられる結果は,できれば python で整理して表示させたい. そのためには,matplotlib を使いこなさなければならない.

しかし,これは python に限ったことではないが,描画機能というのはえてして奥が深い. データを表示させるだけならまだしもそこにはビジュアル要素が絡むから, 多少複雑なことをしようとするとたちどころに覚えなければいけないコマンドが増えたりする. 単純な描画なら plot だけで済むが,人に見せるグラフを描こうと思えば xyラベル追加,凡例追加,範囲指定,ラベルのフォント変更, グラフの種類の変更,mesh の追加,などなど,やることがたくさんある.

幸い python には, ギャラリー という恵まれた例文集が存在する. これを用いつつ,matplotlib の使い方というものを学んでゆきたいと思う.

特異値分解

工学における線形代数に欠かせないのが特異値分解 (singular value decomposition: SVD) である.

特異値分解


任意の行列$A\in \mathbb{C}^{m\times n}$ について, あるユニタリー行列 $U\in \mathbb{C}^{m\times m}, V\in \mathbb{C}^{n\times n}$ が存在し, \begin{align} A = U D V^{\mathsf{H}} \end{align} という形にできる.ただし,$D \in \mathbb{C}^{m\times n}$ は対角成分のみ値を持つ*1


$A$ が正方行列でなおかつエルミート行列の場合, 上記はユニタリー対角化に一致する.このとき,$U=V$ となる. すなわち特異値分解は,エルミート性どころか,正方行列とも保証されていない一般の行列に対して,ユニタリー対角化を拡張したものであるとみなせる.

特異値

固有値分解において, 真ん中の対角行列 $D$ の対角成分が固有値を表すことは良く知られている. 正方行列では,非ゼロ固有値の個数が rank を表すこともまた,線形代数における基礎知識の一つである.

それに対して,特異値分解でできた $D$ の対角成分を特異値と呼ぶ. 非ゼロ特異値の個数は,元の行列 $A$ の rank に等しいことも知られている.

たとえば,$3\times 4$行列 $A$ を特異値分解した結果,[5, -3, 0] という特異値が得られたとする. これは,元の行列の rank が 2 であることを表す.

たとえば,$3\times 4$行列 $A$ を特異値分解した結果,[5, -3, 0.0000001] という特異値が得られたとする. これは,元の行列の rank が 3 であることを表す ... 果たしてそうだろうか.

[5, -3, 0.0000001]というのは確かにどれも非ゼロであるが, 最後の特異値はほかの2つの値に比べて著しく絶対値が小さい. 工学的な話をすると,数値計算の結果出てくる行列では, 数値誤差などにより「厳密に特異」なものはほとんど出てこない. 工学的に特異な行列とは,「特異値が厳密に0」なのではなく, 「特異値がある値 $\epsilon$ より小さい」を表すのである.

特異値分解工学的な rank を知るのに大いに役立つ,といえるだろう.

*1:$D$ は対角行列である,と言いたいが,対角行列は正方行列に定義されるものであるから,ここではそう呼ばないこととする.

Markdown の気づき

記事を書く際には Markdown を使用している.
Markdown の基礎的なことはググればいくらでも出てくるが,
基本的な文法は把握しているという人にとっても,
以下の記事はいろいろな気付きを得ることができると思う.

nanoappli.com

忘れがちな numpy

python で行列演算を行う上で numpy は欠かせない. しかし,numpy.ndarray は行列というよりあくまで配列なので, MATLAB などを使っていた人にとっては奇妙と思える動作をすることがある.

  1. index が 0 から始まる
  2. 1行のみ取り出すと自動的に次元が減る
  3. numpy.ndarray 同士の積や商が elementwise な積や商になる

1. や 3. については,プログラミング言語の配列をイメージすればそこまで不自然なことではない. むしろ,慣れてくれば MATLAB の仕様よりも便利かもしれない. なお,numpy で行列積を計算したければ,@ 演算子を用いるのが便利である.

import numpy as np

A = np.array([[1, 2], [3, 4], [5, 6]])
B = np.array([[2, 3], [4, 5]])

print(A @ B)
  1. についてはどうしてこういう仕様になったのか知らないが,たとえば
A[:,1]

とすると

array([2, 4, 6])

が返ってくる.この shape は (3,) である*1. 本来は(1,3)の形をした縦ベクトルを取り出したいにもかかわらず, (3,)という1次元の配列が返ってきている. なおこの事象は横ベクトルを取り出すときにも同様に起きる.(3,1)を期待しても(3,)の配列が返ってくる.

これを回避したい,すなわち縦ベクトルは縦ベクトルで取り出したい,というときには,大別して2つのやり方がある.

  1. numpy.newaxis または None を使って次元を追加する
  2. スライスを用いる
A[:, np.newaxis, 1]
A[:, [1]]

どちらでも(1,3)の配列が取り出せる.お好みでよさそう.

*1:tuple であることを示すために3の後ろにカンマがついている

このブログについて

数学系,プログラミング系について思いついたことや見知ったことをまとめるブログ.
3行もしないような短い記事を書くかもしれないし,
筆がのったらたくさん書くかもしれない.

プログラミングは,最近 python をやっているので,おそらく python 関係になる.