numpy(「ナンパイ」と読む)は、Pythonで科学技術計算のための中核となるライブラリです。これは、高性能の多次元配列オブジェクトを実装し、高速でメモリ効率のよい処理ツールを提供しています。MATLABに精通している人なら、 Matlabユーザーのための NumPyyという文書から有用な情報が得られるでしょう。
Numpyを使うには、まず次のようにしてnumpyパッケージをimport(取り込む)必要があります:
import numpy as np
まずndarrayの算術演算を紹介します。ndarrayに対する算術操作は、その配列要素すべてに作用します。
data = np.array([[1,2,3],[-1,-2,-3]])
print(data)
print(data * 10)
print(data / 5.0)
print(data + 0.5)
print(data + 0.1*data)
print(data * data)
print(1.0/ data)
統計用の関数も用意されています。meanは平均、varは分散、stdは標準偏差を計算します。これらの用語については「確率と確率分布」で述べます。
# 平均
print(data[0].mean()) # ベクトルに対し計算。 np.mean(data[0]) でもよい
print(data.mean()) # 配列に対し計算。 np.mean(data) でもよい
# 分散
print(data[0].var()) # ベクトルに対し計算。 np.var(data[0]) でもよい
print(data.var()) # 配列に対し計算。 np.var(data) でもよい
# 標準偏差
print(data[0].std()) # ベクトルに対し計算。 np.std(data[0]) でもよい
print(data.std()) # 配列に対し計算。 np.std(data) でもよい
ndarray配列には
print(data.shape) # 2次元配列
print(data.dtype) # 要素は整数
data2 = np.array([1.0, 2.0, 3.0])
print(data2.shape) # 1次元配列
print(data2.dtype) # 要素は浮動小数点数
lst = [1, 2.3, 4.0, 5, 6] # リスト
arr1 = np.array(lst)
print(arr1, arr1.shape, arr1.dtype)
入れ子になったリストからもndarray配列が作れます。
lst2 = [[0,1,2], [3,4,5]]
arr2 = np.array(lst2)
print(arr2, arr2.shape, arr2.dtype)
lst3 = [(0,1,2), (3,4,5)]
arr3 = np.array(lst3)
print(arr3, arr3.shape, arr3.dtype)
lst4 = [(0,1,2), [3,4,5]]
arr4 = np.array(lst4)
print(arr4, arr4.shape, arr4.dtype)
Numpyには、配列を作成するための関数が
関数 | 説明 |
---|---|
array | 引数にリスト、タプル、配列など列挙型のデータをとり、ndarrayを生成する。要素の型は推測、もしくは明示的に指定されたもの |
asarray | arrayと同じくndarrayを生成するが、引数がndarrayの場合は新規には生成しない |
arange | Python組み込みの関数rangeと同じ動作でndarrayを生成する |
ones, ones_like | onesは指定されたサイズのndarrayを、指定されたdtypeで生成し、要素をすべて1とする。ones_likeは引数に別の列挙型オブジェクトを取り、それをテンプレートとして要素をすべて1としたndarray配列を生成する。 |
zeros, zeros_llike | 要素を1ではなく0とすること以外は、onesおよびones_likeと同じ |
empty, empty_like | 要素を初期化せず不定であること以外は、onesおよびones_likeと同じ |
eye, identity | NxNの単位行列となるndarray配列を生成する。 |
a = np.zeros((2,2)) # すべての要素が0である配列を作る
print( a )
b = np.ones((1,2)) # すべての要素が1である配列を作る
print (b)
c = np.full((2,2), 7) # 定数配列を作る
print (c)
d = np.eye(2) # 2x2 の単位行列を作る
print (d)
e = np.random.random((2,2)) # 配列を作り、乱数で要素を埋める
print (e)
配列の別な作成方法については このドキュメントを参照してください。
なお、乱数生成には numpyのrandomモジュールの方が randomモジュールよりもいろいろな関数が用意されています。import random
print(len(dir(random)))
import numpy as np
print(len(dir(np.random)))
arr = np.arange(10) # array([0,1,2,3,4,5,6,7,8,9])
print(arr[5])
ndarrayから切り出した一部(スライス)を取り出すことができます。
print(arr[5:8])
次はarr[5], arr[6], arr[7]の3つの要素それぞれに12という値を代入するものです。これをブロードキャストと呼びます。
arr[5:8] = 12
print(arr)
Pythonのリストとndarrayの違いのひとつは、ndarrayのスライスは元のデータのコピーではなく元のデータを指している、ということです。つまりndarray配列のスライスに対する変更は元のndarray配列のに対する変更になります。これをリストと比較する形で示しましょう。
# Pythonのリスト
lst = list(range(10))
lst_slice = lst[5:8]
lst_slice[1] = 12
print(lst)
# ndarray
arr = np.array(np.arange(10))
arr_slice = arr[5:8]
arr_slice[1] = 12
print(arr)
arr_slice[:]=100
print(arr)
参考: コピーを作るには 次のようにcopy 関数を用います
arr[5:8].copy()
2次元以上の多次元配列の場合、インデックス参照されるものがスカラー値ではなく配列になります。そのためインデックス参照の方法が増えます。
arr2 = np.array([[0,1,2], [3,4,5], [6,7,8]] )
print(arr2[1])
したがって個々の要素を取り出すには、
print(arr2[1][2], arr2[1,2], arr2[(1,2)])
ndarray配列はスライス記法で部分配列を取り出せますが、配列のそれぞれの次元にスライスを指定することができます。
# 次のような階数2で形状(3,4)の配列を作る
# [[ 1 2 3 4]
# [ 5 6 7 8]
# [ 9 10 11 12]]
a = np.array([[1,2,3,4], [5,6,7,8], [9,10,11,12]])
# スライスを用いて最初の2行、1列目と2列目の要素からなる部分配列を抜き出す
# b はそのような形状(2,2)の配列
# [[2 3]
# [6 7]]
b = a[:2, 1:3]
print (b)
配列の真ん中の行にあるデータにアクセスする2通リの方法を見ます。整数インデックスとスライスを混在させることで低い階数の配列ができます。一方、スライスだけを使うと、元の配列と同じ階数の配列ができます。
row_r1 = a[1, :] # 配列aの2行目に対する階数 1 のスライス
row_r2 = a[1:2, :] # 配列aの2行目に対する階数 2 のスライス
row_r3 = a[[1], :] # 配列aの2行目に対する階数 2のスライス
print (row_r1, row_r1.shape )
print (row_r2, row_r2.shape)
print (row_r3, row_r3.shape)
# 配列の列に対しても同じことが起こります
col_r1 = a[:, 1]
col_r2 = a[:, 1:2]
print (col_r1, col_r1.shape)
print (col_r2, col_r2.shape)
スライスを使用してnumpyの配列にインデックスした結果は、常に元の配列の部分配列になります。これと対照的に、整数インデックスは、別の配列のデータを使用して任意の配列を構築することを可能にします。次に例を示します。
a = np.array([[1,2], [3, 4], [5, 6]])
# 整数インデックスの例
# 返される配列は形状(3, )
print (a[[0, 1, 2], [0, 1, 0]] )
# 上の例は次のものと等価
print ( np.array([a[0, 0], a[1, 1], a[2, 0]]) )
# 整数インデックスを使えば、同じ配列の要素を使い回すことができる
print (a[[0, 0], [1, 1]])
# これは前のものと同じ
print (np.array([a[0, 1], a[0, 1]]))
整数配列のインデックスが便利なのは、行列の各行から1つの要素を選択したり変更したりする場合です。
# 要素を選んで配列を作る
a = np.array([[1,2,3], [4,5,6], [7,8,9], [10, 11, 12]])
print (a)
# インデックスの配列を作る
b = np.array([0, 2, 0, 1])
# bのインデックスを用いてそれぞれの行からひとつずつ要素を選ぶ
print (a[np.arange(4), b] ) # "[ 1 6 7 11]"と出力
# bのインデックスを用いてそれぞれの行からひとつずつ要素の値を変化させる
a[np.arange(4), b] += 10
print (a)
ブール配列インデックスを使用すると、配列の任意の要素を取り出すことができます。このタイプのインデックス参照は、ある条件を満たす配列の要素を選択するために頻繁に使用されます。次に例を示します。
import numpy as np
a = np.array([[1,2], [3, 4], [5, 6]])
bool_idx = (a > 2) # 2よりも大きな要素を見つける;
# これによりaと同じ形状のブール値のNumpy配列が返る
# bool_idx のそれぞれの要素は
# 対応するaの要素が2より大きいかどうかの情報を与える
print (bool_idx)
# ブール配列インデックスbool_idxを用いて階数1の配列を作る
# bool_idxの要素がTrueである要素に対応するaの要素からなる
print (a[bool_idx])
# 以上の事を一行で完結に書くことができる:
print (a[a > 2])
話を簡潔にするため、numpy配列のインデックス参照の詳細は省略しました。もっと知りたい場合は このドキュメントを読みましょう。
x = np.array([1, 2]) # Let numpy choose the datatype
y = np.array([1.0, 2.0]) # Let numpy choose the datatype
z = np.array([1, 2], dtype=np.int64) # Force a particular datatype
print (x.dtype, y.dtype, z.dtype)
numpyのデータ型について詳しく知りたい場合は、 このドキュメントのを読むことをお勧めします。
基本的な数学関数は配列の要素ごとに作用し、暗黙的に拡張された演算子や、numpyモジュールの関数によって利用できます。
x = np.array([[1,2],[3,4]], dtype=np.float64)
y = np.array([[5,6],[7,8]], dtype=np.float64)
# 要素ごとの和;両方とも配列を返す
print (x + y)
print (np.add(x, y))
# 要素ごとの差; どちらも配列を返す
print (x - y)
print (np.subtract(x, y))
# 要素ごとの積; どちらも配列を返す
print (x * y)
print (np.multiply(x, y))
# 要素ごとの商; どちらも配列を返す
# [[ 0.2 0.33333333]
# [ 0.42857143 0.5 ]]
print (x / y)
print (np.divide(x, y))
# 要素ごとの平方根; 配列を返す
# [[ 1. 1.41421356]
# [ 1.73205081 2. ]]
print (np.sqrt(x))
*
はMATLABとは異なり、行列乗算ではなく、要素単位の乗算です。numpyでは、dot
関数を使って、行列の積を計算したり、行列とベクトルの積を求めたり、ベクトルの内積を計算します。dot
は、numpyモジュールの関数としても、配列オブジェクトのインスタンス・メソッドとしても利用できます。
x = np.array([[1,2],[3,4]])
y = np.array([[5,6],[7,8]])
v = np.array([9,10])
w = np.array([11, 12])
# ベクトルの内積計算; どちらも 219 を返す
print (v.dot(w))
print (np.dot(v, w))
# 行列/ ベクトルの積; どちらも階数1 の配列 [29 67] を返す
print (x.dot(v))
print (np.dot(x, v))
# 行列 / 行列の積; どちらも階数 2 の配列を返す
# [[19 22]
# [43 50]]
print (x.dot(y))
print (np.dot(x, y))
Numpyは、配列の計算を実行するための多くの便利な関数を提供しています。sum
は最も有用な関数の一つです。
x = np.array([[1,2],[3,4]])
print (np.sum(x)) # すべての要素の総和を計算; "10"と表示する
print (np.sum(x, axis=0)) # 各列の総和を計算; "[4 6]"と表示する
print (np.sum(x, axis=1)) # 各行の総和を計算; "[3 7]"と表示する
numpyが提供する数学関数の完全なリストは このドキュメントにあります。
配列を使用して数学関数を計算することとは別に、配列内のデータを再形成または操作する必要があることがよくあります。このタイプの操作の最も単純な例は行列を転置することですが、
これには単に配列オブジェクトのT
属性を使えばできます:
print (x)
print (x.T)
print(x.transpose()) # transpose関数を使ってもよい
v = np.array([[1,2,3]])
print ( v )
print (v.T)
もっとも行列クラスを使うと、転置行列だけではなく、逆行列を求めることもできる。
a = np.matrix('1 2 ; 3 4') # 行列クラスのインスタンスを作成
# a=np.matrix([[1, 2], [3, 4]]) でも同じ行列が得られる
print(type(a)) # matrixクラス
print(a)
print(a.T) # 転置
inv_a = a.I
print(inv_a) # 逆行列
print(a*inv_a) # 行列の乗算 --- 単位行列になるはず
# ベクトルvを行列xの各行に加算する
# その結果を行列yとして記憶する
x = np.array([[1,2,3], [4,5,6], [7,8,9], [10, 11, 12]])
v = np.array([1, 0, 1])
y = np.empty_like(x) # xと同じサイズの空の行列を作る
# ループを用いてvを行列xの各行に加える
for i in range(4):
y[i, :] = x[i, :] + v
print( y )
これはこれで動くのですが、行列 x
がかなり大きい場合、
Python の機能を使ってこのような繰り返しをやるととても遅くなります。ここで、ベクトル
v
を行列
x
の各行に足すというのは、
v
のコピーを垂直にたくさん積み上げた行列
vv
を作り、それから行列
x
と行列vv
を要素ごとに足すことに等しいことに注目すると、この操作を次のように実現することができます:
vv = np.tile(v, (4, 1)) # vを4つ積み上げる
print (vv) # 出力は: "[[1 0 1]
# [1 0 1]
# [1 0 1]
# [1 0 1]]"
y = x + vv # xとvvを要素ごとに加算する
print (y)
Numpyのブロードキャストにより、 このような計算を、実際にはv
のコピーを作らずにすませることができます。 このような場合にブロードキャストを考えましょう:
import numpy as np
# ベクトルvを行列xの各行に加え
# その結果を行列yで記憶する
x = np.array([[1,2,3], [4,5,6], [7,8,9], [10, 11, 12]])
v = np.array([1, 0, 1])
y = x + v # ブロードキャストを使ってvをxのそれぞれの行に足す
print (y)
y = x + v
は、 たとえx
が
(4, 3)
という形状であり、 v
が(3,)
という形状であっても、ブロードキャストのおかげで計算できます。
あたかもv
が実際に(4, 3)
という
形状をもっていて、各行がみな
v
のコピーであるかのように、和が要素ごとに計算されるのです。
2つの配列に対するブロードキャストは次の規則に従います:
もしもこの説明で意味がわからなければ、 この解説か、 こちらの解説を読んでみよう。
ブロードキャストができる関数は 汎化関数と呼ばれている。どんなものがそうかは、 このドキュメント を参照すること.
以下に、ブロードキャストの使用例をいくつか示す:
import numpy as np
# ベクトルの外積計算
v = np.array([1,2,3]) # v の形状 (3,)
w = np.array([4,5]) # w の形状 (2,)
# 外積を計算するには、まずvの形状を変えて、
# 形状 (3, 1)の列ベクトルにする; 次にこれをwに対してブロードキャストし、
# 形状 (3, 2)の出力を得る。これが v と wの外積である
# [[ 4 5]
# [ 8 10]
# [12 15]]
print (np.reshape(v, (3, 1)) * w)
# ベクトルを行列のそれぞれの「列」に加える
x = np.array([[1,2,3], [4,5,6]])
# x は形状 (2, 3) で v は形状 (3,) であるから、それらをブロードキャストし
# 次のような形状 (2, 3)の行列を得る
# [[2 4 6]
# [5 7 9]]
print (x + v)
# ベクトルを行列の各「行」に加える
# 行列x は形状 (2, 3) で、 w は形状 (2,)である。
# 行列 x を転置すると形状(3, 2) となり、これをwに対しブロードキャストすることで
# 形状(3, 2)の行列を得る;この結果を転置することで
# 形状 (2, 3)の行列を最終結果として得る。これは行列xの各列に
# ベクトル w を加えたものになっており、結果は:
# [[ 5 6 7]
# [ 9 10 11]]
print ( (x.T + w).T )
# 別解: wの形状を形状 (2, 1)のベクトルに変えると
# これを行列 x にブロードキャストして、同じ結果を得る
print (x + np.reshape(w, (2, 1)))
# 行列の定数倍:
# x は形状(2, 3) の行列とする. Numpy ではスカラーは形状 ()のベクトル扱いである
# これらは一緒にブロードキャストされ、形状(2, 3)の次の行列を得る:
# [[ 2 4 6]
# [ 8 10 12]]
print (x * 2)
ブロードキャストを使うと、一般にコードが簡潔でしかも処理速度が速くなります。 ですから、可能な限りこれを使うよう、努めてください。
この短い文書では、知らなければならない重要なことのホンのさわりしか紹介できていないので、 numpyの文献 を読んで、Numpyについて、よりよく知る努力をしてください。
from scipy.misc import imread, imsave, imresize
# JPEG 画像ファイルを読み込み、numpy の配列として記憶する
img = imread('assets/cat.jpg')
print( img.dtype, img.shape ) # "uint8 (400, 248, 3)" と出力する
# 色のチャネルごとにいろいろな定数を掛けてスケーリングすることで
# 画像に色合いをつけることができる。この画像の形状は(400, 248, 3)である
# これに形状(3,) の定数ベクトル[1, 0.95, 0.9] をNumpyのブロードキャスト
# する。つまり赤色(第1チャネル)は変化させないが、緑(第2チャネル)と青
# (第3チャネル)の大きさをそれぞれ0.95倍、および0.9倍する
img_tinted = img * [1, 0.95, 0.9]
# 結果の画像を 300 x 300 ピクセルになるようサイズを変更する
img_tinted = imresize(img_tinted, (300, 300))
# この画像をディスクに書き込む
imsave('assets/cat_tinted.jpg', img_tinted)
scipy.io.loadmat
とscipy.io.savemat
関数を使うと、 MATLABのファイルの読み書きができます。
これについては
このドキュメントをご覧ください.
SciPyは、点(座標)の集合の間の距離を計算するための有用な関数を用意しています、
関数 scipy.spatial.distance.pdist
は、
集合中のすべての点と点のペアの間の距離を計算します
import numpy as np
from scipy.spatial.distance import pdist, squareform
# 各行が2D空間における点となっている配列を作る:
# [[0 1]
# [1 0]
# [2 0]]
x = np.array([[0, 1], [1, 0], [2, 0]])
print (x)
# xのすべての行の間のユークリッド距離を計算する
# d[i, j] は x[i, :] と x[j, :] の間の距離となる
# 結果は次:
# [[ 0. 1.41421356 2.23606798]
# [ 1.41421356 0. 1. ]
# [ 2.23606798 1. 0. ]]
d = squareform(pdist(x, 'euclidean'))
print (d)
Matplotlibは、描画のためのライブラリです。
ここでは、 matplotlib.pyplot
モジュールについて簡単に紹介します。
これは、MATLABの描画機能と類似のシステムを提供しています。
matplotlibで最も重要な関数が plot
です。これにより
2Dデータをプロットすることができます。例を示します:
%matplotlib inline
# 上により、matplotlibで作成した図がブラウザ上に表示されるようになる
import numpy as np
import matplotlib.pyplot as plt
# x座標に対するsin関数の値を計算しy座標の値とする
x = np.arange(0, 3 * np.pi, 0.05)
y = np.sin(x)
# x軸を引く
plt.axhline(y=0, xmin=min(x), xmax=max(x), color='green',linestyle='dashed')
# matplotlibを使ってプロットする
plt.plot(x, y)
plt.show() # 本来はこの plt.show() を呼ばなければ図が見えない
# sin関数を散布図として表示する
x = np.arange(0, 3 * np.pi, 0.05)
y = np.sin(x)
plt.plot(x,y,'b.') # 第3引数に色(Blue)とマーカーの種類を指定すると散布図となる
plt.grid() # 枠線を書く
plt.show()
ちょっとだけ余計に作業することで、簡単にいくつもの線を同時に描画したり、題目をつけたり、 凡例をつけたり、軸のラベルをつけたりできます:
import numpy as np
import matplotlib.pyplot as plt
# x座標に対するsin関数の値を計算しy座標の値とする
x = np.arange(0, 3 * np.pi, 0.1)
y_sin = np.sin(x)
y_cos = np.cos(x)
# matplotlibを使ってプロットする
plt.plot(x, y_sin)
plt.plot(x, y_cos)
plt.xlabel('x axis label') # ラベルに日本語は使えないと思って良い
plt.ylabel('y axis label')
plt.title('Sine and Cosine')
plt.legend(['Sine', 'Cosine'])
plt.grid()
plt.show() # これを最後に呼ぶこと
plot
関数について、詳しくは
この文書をご覧ください
import numpy as np
import matplotlib.pyplot as plt
# x座標に対応するsinとcos関数の値を計算する
x = np.arange(0, 3 * np.pi, 0.1)
y_sin = np.sin(x)
y_cos = np.cos(x)
# subplotの枠を設定:高さ2、幅1とする,
# そして、最初の枠をactiveにする
plt.subplot(2, 1, 1)
# 最初の図を描画する
plt.plot(x, y_sin)
plt.title('Sine')
plt.grid()
# 2番めの枠をactiveにして、別な図を描画する
plt.subplot(2, 1, 2)
plt.plot(x, y_cos)
plt.title('Cosine')
plt.grid()
# 図を表示する(これを最後に呼ぶこと)
plt.show()
subplot
関数について、詳しくは
この文書をご覧ください
imshow
関数を使って画像を表示できます。次の例をみてください。ただし以下を実行するには
assets.zipをダウンロードし展開しておくこと
import numpy as np
from scipy.misc import imread, imresize
import matplotlib.pyplot as plt
img = imread('assets/cat.jpg')
img_tinted = img * [1, 0.95, 0.9]
# 原画像を表示
plt.subplot(1, 2, 1)
plt.imshow(img)
# 色を変えた画像を表示
plt.subplot(1, 2, 2)
# A slight gotcha with imshow is that it might give strange results
# if presented with data that is not uint8. To work around this, we
# explicitly cast the image to uint8 before displaying it.
plt.imshow(np.uint8(img_tinted))
plt.show()
# 問題1
# 問題2
3. 行列の固有値と固有ベクトルは次のようにして求めることができる(参考:固有値、固有ベクトルを可視化する)
# Mを行列とする
eigenVals, eigenVecs = np.linalg.eig(M) # eigenValsが固有値、eigenVecsが固有ベクトルの行列
これに基づき、Aの値である行列の固有値と固有ベクトルを求めよ。(注意: 固有ベクトルは、固有ベクトルの行列の「列」ベクトルである)
# 問題3
4. 3x3の行列に対し、固有値とそれに対応する固有ベクトルは三組ある。
元の行列を$\mathbf{A}$、 固有値と対応する固有ベクトルをそれぞれ $\lambda$と$\mathbf{x}$ と書くと、 $\mathbf{Ax}=\lambda\mathbf{x}$ が成り立つこと、また固有ベクトルがみな正則化されている(長さが1)であることも確認せよ。
ただし固有値が実数のものだけに限定するとする(つまり固有値が複素数のものは除外する)・
# 問題4
5. 次を実行し、散布図をmatplotlibを用いて以下のような図を表示せよ。
x=np.array([2.8, 2.9, 3.0, 3.1, 3.2, 3.2, 3.2, 3.3, 3.4]) y=np.array([30,26,33,31,33,35,37,36,33])
# 問題5
6. 最小二乗法を用いて、上記のxとyの関係を近似する線形関数$y = ax + b$を求めよう(つまり、傾き$a$と切片$b$を求める)。
これは、上記のデータを$x_i, y_i$ (i=1,..,9)とおくと、
$ \left( \begin{array}{r} a \\ b \end{array} \right) =
\left( \begin{array}{rr} \sum_i x_i^2 & \sum_i x_i \\ \sum_i x_i & \sum_i 1 \\ \end{array} \right)^{-1}
\left( \begin{array}{r} \sum_i x_i y_i \\ \sum_i y_i \end{array} \right) $
という式から得られる。これによって近似線形関数を求め、上の図にその直線を書き込め。次のような図が得られれば良い。
注意: 上記を満たす関数は後で学ぶ。ここでは上記の方法にしたがって解いてみよ。
直線は両端の2点を与えてplotすればよい。例えばplt.plot([0,1],[0,10],'g')
で(0,0)と(1,10)を結ぶ直線がひける
# 問題6