標準インストールの Python では、リスト (list) を使って一連の値を保存できますが、これは配列として使用することができます。しかし、リストの要素は任意のオブジェクトであるため、リストに保存されるのはオブジェクトのポインタです。このため、単純な [1,2,3] を保存するには、3 つのポインタと 3 つの整数オブジェクトが必要です。数値計算において、この構造は明らかにメモリと CPU 計算時間を無駄にしています。
さらに、Python は array モジュールも提供しており、array オブジェクトはリストとは異なり、数値を直接保存し、C 言語の 1 次元配列に似ています。しかし、多次元をサポートしておらず、さまざまな演算関数もないため、数値計算には適していません。
NumPy の誕生はこれらの欠点を補い、NumPy は 2 つの基本オブジェクトを提供します:ndarray(N 次元配列オブジェクト)と ufunc(ユニバーサル関数オブジェクト)。ndarray(以下、配列と呼ぶ)は単一のデータ型を持つ多次元配列であり、ufunc は配列を処理するための関数です。#
ここで、リスト(list)、タプル(tuple)、および配列(array)の違いについて説明します:
リスト: 要素は任意のタイプのオブジェクトであり、初期化が必要で、初期化後の内容は変更可能です。すべての要素オブジェクトが同じである場合、配列として使用できます。
>>> a = ['a', "hello", 1, 2.2]
>>> a
['a', 'hello', 1, 2.2]
>>> a[1] = "HELLO"
>>> a
['a', 'HELLO', 1, 2.2]
タプル: 要素は任意のタイプのオブジェクトであり、初期化が必要で、初期化後の内容は変更不可です。
>>> b = ('a', "hello", 1, 2.2)
>>> b
('a', 'hello', 1, 2.2)
>>> b[2] = "HELLO"
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
TypeError: 'tuple' object does not support item assignment
配列: 要素は 1 つのタイプのオブジェクトのみであり、初期化は必須ではなく、初期化後の内容は変更可能です。
>>> c = np.array([1,2,3])
>>> d = np.array(['a','b','c'])
>>> e = np.array([1,2,'c','d'])
>>> e
array(['1', '2', 'c', 'd'],
dtype='|S21')
見ての通り、配列 e では 1,2 が文字列型として扱われ、後の‘c'、’d' と一致しています。
観察力があれば、タプルは小括弧()、リストと配列は中括弧 [] を使用していることがわかります。#
1. 関数ライブラリのインポート#
import numpy as np
numpy 関数ライブラリをインポートし、以降は略称 np を使用します。
2. 作成#
a、b は 1 次元配列で、c は多次元配列です。
>>> a = np.array([1, 2, 3, 4])
>>> b = np.array((5, 6, 7, 8))
>>> c = np.array([[1, 2, 3, 4],[4, 5, 6, 7], [7, 8, 9, 10]])
>>> b
array([5, 6, 7, 8])
>>> c
array([[1, 2, 3, 4],
[4, 5, 6, 7],
[7, 8, 9, 10]])
>>> c.dtype
dtype('int32')
配列のサイズはその shape 属性を通じて取得できます:
>>> a.shape
(4,)
>>> c.shape
(3, 4)
配列 a の shape は 1 つの要素しか持たないため、1 次元配列です。一方、配列 c の shape は 2 つの要素を持つため、2 次元配列です。0 軸(3 つの横軸)の長さは 3、1 軸(4 つの縦軸)の長さは 4 です。また、配列の shape 属性を変更することで、配列の要素数を変えずに各軸の長さを変更できます。以下の例では、配列 c の shape を (4,3) に変更します。注意すべきは、(3,4) から (4,3) に変更することは配列を転置するのではなく、各軸のサイズを変更するだけで、配列の要素がメモリ内での位置は変わらないことです:
>>> c.shape = 4,3
>>> c
array([[ 1, 2, 3],
[ 4, 4, 5],
[ 6, 7, 7],
[ 8, 9, 10]])
上記の配列 a は 1 次元配列で、a.shape は 1 つの数値しか持たない。この数値は配列内の要素の数を示します。多次元配列の表現方法を使って a を表現することもできます。違いを観察してください:
>>> a = np.array([1, 2, 3, 4])
>>> a
array([1, 2, 3, 4])
>>> a.shape
(4,)
>>> a.shape = 1,4
>>> a
array([[1, 2, 3, 4]])
>>> a.shape
(1, 4)
ある軸の要素が - 1 の場合、その軸の長さは配列の要素数に基づいて自動的に計算されます。したがって、以下のプログラムは配列 c の shape を (2,6) に変更します:
>>> c.shape = 2,-1
>>> c
array([[ 1, 2, 3, 4, 4, 5],
[ 6, 7, 7, 8, 9, 10]])
配列の reshape メソッドを使用すると、サイズを変更した新しい配列を作成できます。元の配列の shape は変わりません:
>>> d = a.reshape((2,2))
>>> d
array([[1, 2],
[3, 4]])
>>> a
array([1, 2, 3, 4])
配列 a と d は実際にはデータストレージのメモリ領域を共有しているため、いずれかの配列の要素を変更すると、もう一方の配列の内容も同時に変更されます:
>>> a[1] = 100 # 配列aの第1要素を100に変更
>>> d # 配列dの2も変更されたことに注意
array([[ 1, 100],
[ 3, 4]])
配列の要素タイプは dtype 属性を通じて取得できます。上記の例では、パラメータシーケンスの要素はすべて整数であるため、作成された配列の要素タイプも整数であり、32 ビットの長整数型です。dtype パラメータを使用して、作成時に要素タイプを指定できます:
>>> np.array([[1, 2, 3, 4],[4, 5, 6, 7], [7, 8, 9, 10]], dtype=np.float)
array([[ 1., 2., 3., 4.],
[ 4., 5., 6., 7.],
[ 7., 8., 9., 10.]])
>>> np.array([[1, 2, 3, 4],[4, 5, 6, 7], [7, 8, 9, 10]], dtype=np.complex)
array([[ 1.+0.j, 2.+0.j, 3.+0.j, 4.+0.j],
[ 4.+0.j, 5.+0.j, 6.+0.j, 7.+0.j],
[ 7.+0.j, 8.+0.j, 9.+0.j, 10.+0.j]])
上記はすべて array 関数を使用して配列を作成しています。以下は配列を作成するための 3 つの専用関数です。
arange 関数は Python の range 関数に似ており、開始値、終了値、およびステップを指定して 1 次元配列を作成します。 注意:配列には終了値は含まれません :
>>> np.arange(0,1,0.1)
array([ 0. , 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9])
linspace 関数は開始値、終了値、および要素数を指定して 1 次元配列を作成します。endpoint キーワードを使用して終了値を含むかどうかを指定できます。デフォルト設定は 終了値を含む です :
>>> np.linspace(0, 1, 12)
array([ 0. , 0.09090909, 0.18181818, 0.27272727, 0.36363636,
0.45454545, 0.54545455, 0.63636364, 0.72727273, 0.81818182,
0.90909091, 1. ])
logspace 関数は linspace に似ていますが、等比数列を作成します。以下の例は 1 (10^0) から 100 (10^2) までの 20 個の要素を持つ等比数列を生成します:
>>> np.logspace(0, 2, 20)
array([ 1. , 1.27427499, 1.62377674, 2.06913808,
2.6366509 , 3.35981829, 4.2813324 , 5.45559478,
6.95192796, 8.8586679 , 11.28837892, 14.38449888,
18.32980711, 23.35721469, 29.76351442, 37.92690191,
48.32930239, 61.58482111, 78.47599704, 100. ])
さらに、frombuffer、fromstring、fromfile などの関数を使用してバイトシーケンスから配列を作成できます。以下は fromstring の例です:
>>> s = "abcdefgh"
Python の文字列は実際にはバイトシーケンスであり、各文字は 1 バイトを占めるため、文字列 s から 8 ビットの整数配列を作成すると、得られる配列は文字列中の各文字の ASCII コードになります:
>>> np.fromstring(s, dtype=np.int8)
array([ 97, 98, 99, 100, 101, 102, 103, 104], dtype=int8)
もし文字列 s から 16 ビットの整数配列を作成すると、隣接する 2 つのバイトが 1 つの整数を表し、バイト 98 とバイト 97 を 1 つの 16 ビット整数として扱うと、その値は 98*256+97 =
25185 になります。メモリ内ではリトルエンディアン(低位バイトが前)方式でデータが保存されていることがわかります。
>>> np.fromstring(s, dtype=np.int16)
array([25185, 25699, 26213, 26727], dtype=int16)
>>> 98*256+97
25185
Python の関数を作成して、配列のインデックスを配列内の対応する値に変換し、この関数を使用して配列を作成できます:
>>> def func(i): #関数を定義、関数名func、引数i
... return i%4+1 #iを4で割った余りに1を加えた値を返す
...
>>> np.fromfunction(func, (10,))
array([ 1., 2., 3., 4., 1., 2., 3., 4., 1., 2.])
fromfunction の第 2 引数 n は、0 から始まる n 個の値を取得し、前の function に割り当てられます。
この例では、i が最初に 0 を取得し、余りに 1 を加えた後、値 1 を返します。
次に i が 1 を取得し、余りに 1 を加えた後、値 2 を返します。
このように続きます。
最後に i が 9 を取得し、余りに 1 を加えた後、値 2 を返します。#
fromfunction 関数の第 1 引数は各配列要素を計算する関数であり、第 2 引数は配列のサイズ (shape) です。多次元配列をサポートしているため、第 2 引数はシーケンスでなければなりません。この例では (10,) を使用して 10 要素の 1 次元配列を作成します。
以下の例では、九九の掛け算表を表す 2 次元配列を作成し、出力される配列 a の各要素 a [i, j] は func2 (i, j) に等しいです:
>>> def func2(i, j): #関数を定義、関数名func、引数i、j
... return (i+1) * (j+1) #関数の返り値
...
>>> a = np.fromfunction(func2, (9,9))
>>> a
array([[ 1., 2., 3., 4., 5., 6., 7., 8., 9.],
[ 2., 4., 6., 8., 10., 12., 14., 16., 18.],
[ 3., 6., 9., 12., 15., 18., 21., 24., 27.],
[ 4., 8., 12., 16., 20., 24., 28., 32., 36.],
[ 5., 10., 15., 20., 25., 30., 35., 40., 45.],
[ 6., 12., 18., 24., 30., 36., 42., 48., 54.],
[ 7., 14., 21., 28., 35., 42., 49., 56., 63.],
[ 8., 16., 24., 32., 40., 48., 56., 64., 72.],
[ 9., 18., 27., 36., 45., 54., 63., 72., 81.]]) #各点のインデックスを掛け算したもの
3. 要素の取得と保存#
配列要素の取得方法は Python の標準的な方法と同じです:
>>> a = np.arange(10)
>>> a[5] # 整数をインデックスとして使用して配列内の特定の要素を取得
5
>>> a[3:5] # 範囲をインデックスとして使用して配列のスライスを取得、a[3]を含みa[5]を含まない
array([3, 4])
>>> a[:5] # 開始インデックスを省略し、a[0]から始めることを示す
array([0, 1, 2, 3, 4])
>>> a[:-1] # インデックスに負の数を使用し、配列の最後の要素から前に数える
array([0, 1, 2, 3, 4, 5, 6, 7, 8])
>>> a[2:4] = 100,101 # インデックスを使用して要素の値を変更
>>> a
array([ 0, 1, 100, 101, 4, 5, 6, 7, 8, 9])
>>> a[1:-1:2] # 範囲の第3引数はステップを示し、2は1つおきに要素を取得することを示す
array([ 1, 101, 5, 7])
>>> a[::-1] # 範囲の開始インデックスと終了インデックスを省略し、ステップを-1に設定して配列全体を逆にする
array([ 9, 8, 7, 6, 5, 4, 101, 100, 1, 0])
>>> a[5:1:-2] # ステップが負の数の場合、開始インデックスは終了インデックスより大きくなければならない
array([ 5, 101])
Python のリストシーケンスとは異なり、インデックス範囲を使用して取得した新しい配列は元の配列のビューです。元の配列と同じデータ空間を共有しています:
>>> b = a[3:7] # インデックス範囲を使用して新しい配列bを生成し、bとaは同じデータ空間を共有
>>> b
array([101, 4, 5, 6])
>>> b[2] = -10 # bの第2要素を-10に変更
>>> b
array([101, 4, -10, 6])
>>> a # aの第5要素も変更されている
array([ 0, 1, 100, 101, 4, -10, 6, 7, 8, 9])
インデックス範囲を使用して要素を取得する以外に、NumPy は要素を取得するための 2 つの高度な方法も提供しています。
整数シーケンスを使用
整数シーケンスを使用して配列要素を取得する場合、シーケンス内の各要素をインデックスとして使用します。整数シーケンスはリストまたは配列であることができます。整数シーケンスをインデックスとして使用して取得した配列は元の配列とデータ空間を共有しません。
>>> x = np.arange(10,1,-1)
>>> x
array([10, 9, 8, 7, 6, 5, 4, 3, 2])
>>> x[[3, 3, 1, 8]] # xのインデックス3, 3, 1, 8の4つの要素を取得し、新しい配列を構成
array([7, 7, 9, 2])
>>> b = x[np.array([3,3,-3,8])] #インデックスは負の数でも構いません
>>> b[2] = 100
>>> b
array([7, 7, 100, 2])
>>> x # bとxはデータ空間を共有しないため、xの値は変更されていません
array([10, 9, 8, 7, 6, 5, 4, 3, 2])
>>> x[[3,5,1]] = -1, -2, -3 # 整数シーケンスインデックスも要素の値を変更するために使用できます
>>> x
array([10, -3, 8, -1, 6, -2, 4, 3, 2])
ブール配列を使用
ブール配列 b をインデックスとして使用して配列 x の要素を取得する場合、配列 b で対応するインデックスが True である配列 x のすべての要素を収集します。ブール配列をインデックスとして使用して取得した配列は元の配列とデータ空間を共有しません。この方法はブール配列にのみ対応し、ブールリストを使用することはできません。(リスト、タプル、配列の違いを覚えていますか?)
>>> x = np.arange(5,0,-1)
>>> x
array([5, 4, 3, 2, 1])
>>> x[np.array([True, False, True, False, False])]
>>> # ブール配列のインデックス0,2の要素がTrueであるため、xのインデックス0,2の要素を取得
array([5, 3])
>>> x[[True, False, True, False, False]]
>>> # ブールリストの場合、Trueを1、Falseを0として扱い、整数シーケンス方式でxの要素を取得
array([4, 5, 4, 5, 5])
>>> x[np.array([True, False, True, True])]
>>> # ブール配列の長さが不足している場合、不足分はすべてFalseとして扱われる
array([5, 3, 2])
>>> x[np.array([True, False, True, True])] = -1, -2, -3
>>> # ブール配列インデックスも要素を変更するために使用できます
>>> x
array([-1, 4, -2, -3, 1])
ブール配列は通常手動で生成されるのではなく、ブール演算の ufunc 関数を使用して生成されます。ufunc 関数については「ufunc 運算」セクションを参照してください。
>>> x = np.random.rand(10) # 長さ10、要素値が0-1のランダム数の配列を生成
>>> x
array([ 0.72223939, 0.921226 , 0.7770805 , 0.2055047 , 0.17567449,
0.95799412, 0.12015178, 0.7627083 , 0.43260184, 0.91379859])
>>> x>0.5
>>> # 配列xの各要素と0.5を比較し、ブール配列を生成。Trueはxの対応する値が0.5より大きいことを示す
array([ True, True, True, False, False, True, False, True, False, True], dtype=bool)
>>> x[x>0.5]
>>> # x>0.5から返されたブール配列を使用してxの要素を収集し、結果として得られるのはxのすべての要素の配列
array([ 0.72223939, 0.921226 , 0.7770805 , 0.95799412, 0.7627083 ,
0.91379859])
4. 多次元配列#
多次元配列の取得方法は 1 次元配列と似ています。多次元配列には複数の軸があるため、そのインデックスは複数の値で表す必要があります。NumPy はタプル (tuple) を配列のインデックスとして使用します。以下の図では、a は 6x6 の配列で、図中の色で各インデックスとその対応する選択領域を区別しています。
a [2::2,::2] という文を分解して見ると、最初の「2:」はインデックス 2 の横軸(0 軸)から末尾までを示し、「:2」は 2 行おきに取得することを示します。つまり、2 行目と 4 行目を取得します。後の最初の「:」はすべての列(1 軸)を示し、「:2」は 2 列おきに取得することを示します。つまり、0 列、2 列、4 列を取得します。
例えば 24 のインデックスは (2,4) で、2 は 0 軸のインデックス、4 は 1 軸のインデックスです。
配列 a は加算表で、縦軸の値は 0, 10, 20, 30, 40, 50;横軸の値は 0, 1, 2, 3, 4, 5 です。縦軸の各要素は横軸の各要素と加算されます。
>>> np.arange(0, 60, 10).reshape(-1, 1) + np.arange(0, 6)
array([[ 0, 1, 2, 3, 4, 5],
[10, 11, 12, 13, 14, 15],
[20, 21, 22, 23, 24, 25],
[30, 31, 32, 33, 34, 35],
[40, 41, 42, 43, 44, 45],
[50, 51, 52, 53, 54, 55]])
多次元配列も整数シーケンスやブール配列を使用して取得できます。
a[(0,1,2,3,4),(1,2,3,4,5)] :
配列のインデックスを取得するために、2 つの要素を持つタプルを使用し、タプル内の各要素は整数シーケンスであり、それぞれ配列の第 0 軸と第 1 軸に対応します。2 つのシーケンスの対応する位置から 2 つの整数を取得してインデックスを構成します:
a[0,1], a[1,2], ..., a[4,5]。
a [3:, [0, 2, 5]] : インデックスの第 0 軸は範囲であり、第 3 行以降のすべての行を選択します;第 1 軸は整数シーケンスであり、0, 2, 5 の 3 列を選択します。
a [mask, 2] : インデックスの第 0 軸はブール配列であり、0, 2, 5 行を選択します;第 1 軸は整数であり、2 列を選択します。
5. 配列構造#
構造配列を定義する必要があると仮定します。各要素には name、age、weight フィールドがあります。NumPy では次のように定義できます:
import numpy as np
persontype = np.dtype({
'names':['name', 'age', 'weight'],
'formats':['S32','i', 'f']})
a = np.array([("Zhang",32,75.5),("Wang",24,65.2)],
dtype=persontype)
まず、dtype オブジェクト persontype を作成し、その辞書パラメータで構造タイプの各フィールドを説明します。辞書には 2 つのキーワードがあります:names、formats。各キーワードに対応する値はリストです。names は構造内の各フィールド名を定義し、formats は各フィールドのタイプを定義します:
S32 : 32 バイトの文字列型で、構造内の各要素のサイズは固定である必要があるため、文字列の長さを指定する必要があります。
i : 32 ビットの整数型で、np.int32 に相当します。
f : 32 ビットの単精度浮動小数点数型で、np.float32 に相当します。
次に、array 関数を呼び出して配列を作成し、キーワードパラメータ dtype=persontype を指定して、作成された配列の要素タイプを構造 persontype に指定します。上記のプログラムを実行した後、IPython で次のような文を実行して配列 a の要素タイプを確認できます。
>>> a.dtype
dtype([('name', '|S32'), ('age', '<i4'), ('weight', '<f4')])
構造配列の取得方法は一般の配列と同じで、インデックスを使用してその要素を取得できます。要素の値はタプルのように見えますが、実際には構造です:
>>> a[0]
('Zhang', 32, 75.5)
>>> a[0].dtype
dtype([('name', '|S32'), ('age', '<i4'), ('weight', '<f4')])
a [0] は構造要素であり、配列 a とメモリデータを共有しているため、そのフィールドを変更することで元の配列の対応するフィールドを変更できます:
>>> c = a[1]
>>> c["name"] = "Li"
>>> a[1]["name"]
"Li"
構造は辞書のように文字列インデックスを使用して対応するフィールド値を取得できます:
>>> a[0]["name"]
'Zhang'
構造要素の特定のフィールドを取得するだけでなく、構造配列のフィールドを直接取得することもでき、これは元の配列のビューを返します。したがって、b [0] を変更することで a [0]["age"] を変更できます:
>>> b=a[:]["age"] # またはa["age"]
>>> b
array([32, 24])
>>> b[0] = 40
>>> a[0]["age"]
40
a.tostring または a.tofile メソッドを呼び出すことで、配列 a のバイナリ形式を直接出力できます:
>>> a.tofile("test.bin")
6. ufunc 運算#
ufunc はユニバーサル関数の略で、配列の各要素に操作を行う関数です。NumPy に組み込まれている多くの ufunc 関数は C 言語レベルで実装されているため、計算速度が非常に速いです。例を見てみましょう:
>>> x = np.linspace(0, 2*np.pi, 10)
# 配列xの各要素に対して正弦計算を行い、同じサイズの新しい配列を返します
>>> y = np.sin(x)
>>> y
array([ 0.00000000e+00, 6.42787610e-01, 9.84807753e-01,
8.66025404e-01, 3.42020143e-01, -3.42020143e-01,
-8.66025404e-01, -9.84807753e-01, -6.42787610e-01,
-2.44921271e-16])
最初に linspace を使用して 0 から 2*PI までの等間隔の 10 個の数を生成し、それを sin 関数に渡します。np.sin は ufunc 関数であるため、x の各要素に対して正弦値を計算し、その結果を返し、y に代入します。計算後、x の値は変更されず、新しい配列が結果を保存します。もし sin 関数が計算した結果を直接配列 x に上書きしたい場合は、上書きされる配列を ufunc 関数の第 2 引数として渡すことができます。例えば:
>>> t = np.sin(x,x)
>>> x
array([ 0.00000000e+00, 6.42787610e-01, 9.84807753e-01,
8.66025404e-01, 3.42020143e-01, -3.42020143e-01,
-8.66025404e-01, -9.84807753e-01, -6.42787610e-01,
-2.44921271e-16])
>>> id(t) == id(x)
True
sin 関数の第 2 引数も x であるため、x の各値に対して正弦値を計算し、その結果を x の対応する位置に保存します。この時、関数の返り値は計算全体の結果であり、x そのものであるため、2 つの変数の id は同じです(変数 t と変数 x は同じメモリ領域を指しています)。
以下の小さなプログラムでは、numpy.math と Python 標準ライブラリの math.sin の計算速度を比較しています:
import time
import math
import numpy as np
x = [i * 0.001 for i in xrange(1000000)]
start = time.clock()
for i, t in enumerate(x):
x[i] = math.sin(t)
print "math.sin:", time.clock() - start
x = [i * 0.001 for i in xrange(1000000)]
x = np.array(x)
start = time.clock()
np.sin(x,x)
print "numpy.sin:", time.clock() - start
# 出力
# math.sin: 1.15426932753
# numpy.sin: 0.0882399858083
100 万回の正弦値計算を行うと、numpy.sin は math.sin よりも 10 倍以上速いです。これは numpy.sin が C 言語レベルでのループ計算を行うためです。numpy.sin は単一の数値に対しても正弦を計算することができますが、注意すべきは、単一の数値の計算では math.sin の方が速いということです。以下のテストプログラムを見てみましょう:
x = [i * 0.001 for i in xrange(1000000)]
start = time.clock()
for i, t in enumerate(x):
x[i] = np.sin(t)
print "numpy.sin loop:", time.clock() - start
# 出力
# numpy.sin loop: 5.72166965355
numpy.sin の計算速度は math.sin の 1/5 に過ぎません。これは numpy.sin が配列と単一の値の計算を同時にサポートするため、その C 言語の内部実装が math.sin よりも複雑であるためです。もし同様に Python レベルでループを行った場合、その違いがわかります。また、numpy.sin が返す数のタイプは math.sin が返すタイプとは異なります。math.sin は Python の標準 float 型を返しますが、numpy.sin は numpy.float64 型を返します:
>>> type(math.sin(0.5))
<type 'float'>
>>> type(np.sin(0.5))
<type 'numpy.float64'>
上記の例を通じて、math ライブラリと numpy ライブラリの数学関数を最も効率的に使用する方法を理解しました。
NumPy にはさまざまな ufunc 関数があり、さまざまな計算を提供しています。単一入力関数の sin の他にも、多数の入力を持つ関数があり、add 関数は最も一般的な例です。まずは例を見てみましょう:
>>> a = np.arange(0,4)
>>> a
array([0, 1, 2, 3])
>>> b = np.arange(1,5)
>>> b
array([1, 2, 3, 4])
>>> np.add(a,b)
array([1, 3, 5, 7])
>>> np.add(a,b,a)
array([1, 3, 5, 7])
>>> a
array([1, 3, 5, 7])
add 関数は新しい配列を返し、この配列の各要素は 2 つのパラメータ配列の対応する要素の和です。第 3 引数を指定すると、計算結果をどの配列に書き込むかを指定できます。指定した場合、add 関数は新しい配列を生成しません。
Python の演算子オーバーロード機能により、2 つの配列の加算は単純に a+b と書くことができ、np.add (a,b,a) は a+=b と表すことができます。以下は配列の演算子とそれに対応する ufunc 関数のリストです。除算 "/" の意味は__future__.division が有効かどうかによって異なることに注意してください。
y = x1 + x2: add(x1, x2 [, y])
y = x1 - x2: subtract(x1, x2 [, y])
y = x1 * x2: multiply (x1, x2 [, y])
y = x1 /x2: divide (x1, x2 [, y]), もし 2 つの配列の要素が整数の場合、整数除法を使用
y = x1 /x2: true divide (x1, x2 [, y]), 常に正確な商を返す
y = x1 //x2: floor divide (x1, x2 [, y]), 常に返り値を切り捨てる
y = -x: negative(x [,y])
y = x1**x2: power(x1, x2 [, y])
y = x1 % x2: remainder(x1, x2 [, y]), mod(x1, x2, [, y])
配列オブジェクトはこれらの操作をサポートしており、式の記述を大幅に簡素化します。ただし、式が非常に複雑で、計算する配列が非常に大きい場合、多くの中間結果が生成されるため、プログラムの計算効率が低下する可能性があります。例えば、a、b、c の 3 つの配列が式 x=a*b+c を使用して計算される場合、これは次のようになります:
t = a * b
x = t + c
del t
つまり、乗算の計算結果を保存するために配列 t を生成する必要があり、その後最終結果配列 x を生成します。式を手動で分解して x = a*b; x += c とすることで、メモリ割り当てを 1 回減らすことができます。
標準の ufunc 関数の呼び出しを組み合わせることで、さまざまな式の配列計算を実現できます。しかし、時にはこのような式を書くのが難しい場合があり、各要素の計算関数を Python で簡単に実装できる場合があります。このとき、frompyfunc 関数を使用して単一要素の計算関数を ufunc 関数に変換できます。これにより、生成された ufunc 関数を使用して配列を計算することができます。例を見てみましょう。
三角波を記述するための分段関数を使用したいと考えています。三角波の形は以下のようになります:
上記の図に基づいて、三角波のある点の y 座標を計算する関数を次のように書くことができます:
def triangle_wave(x, c, c0, hc):
x = x - int(x) # 三角波の周期は1であるため、x座標の小数部分のみを計算に使用
if x >= c: r = 0.0
elif x < c0: r = x / c0 * hc
else: r = (c-x) / (c-c0) * hc
return r
明らかに triangle_wave 関数は単一の数値しか計算できず、配列を直接処理することはできません。
別の方法もあります:
def triangle_func(c, c0, hc):
def trifunc(x):
x = x - int(x) # 三角波の周期は1であるため、x座標の小数部分のみを計算に使用
if x >= c: r = 0.0
elif x < c0: r = x / c0 * hc
else: r = (c-x) / (c-c0) * hc
return r
# trifunc関数を使用してufunc関数を作成し、配列に対して直接計算できます。ただし、この関数で計算した結果は
# Object配列となるため、型変換が必要です
return np.frompyfunc(trifunc, 1, 1)
y2 = triangle_func(0.6, 0.4, 1.0)(x)
triangle_func 内で三角波の 3 つのパラメータをラップし、その内部で三角波を計算する関数 trifunc を定義します。trifunc 関数は呼び出すと triangle_func のパラメータを使用して計算を行います。最後に triangle_func は frompyfunc を使用して結果を返します。
frompyfunc で得られた関数で計算された配列要素のタイプは object であることに注意してください。これは frompyfunc 関数が Python 関数が返すデータ型が完全に一致することを保証できないためです。したがって、再度 y2.astype (np.float64) を使用してそれをダブル精度浮動小数点配列に変換する必要があります。
7. ブロードキャスト#
ufunc 関数を使用して 2 つの配列を計算する際、ufunc 関数はこれら 2 つの配列の対応する要素を計算します。したがって、これら 2 つの配列は同じサイズである必要があります(shape が同じである必要があります)。もし 2 つの配列の shape が異なる場合、次のようにブロードキャスト (broadcasting) 処理が行われます:
1. すべての入力配列は shape が最も長い配列に合わせられ、shape が不足している部分は前に 1 を追加して補完されます。
2. 出力配列の shape は入力配列の shape の各軸の最大値です。
3. 入力配列のある軸の長さが出力配列の対応する軸の長さと同じか、長さが 1 の場合、その配列は計算に使用できます。そうでない場合はエラーになります。
4. 入力配列のある軸の長さが 1 の場合、その軸に沿って計算する際は、その軸の最初の値を使用します。
上記の 4 つのルールは理解するのが難しいかもしれませんが、実際の例を見てみましょう。
まず、shape が (6,1) の 2 次元配列 a を作成します:
>>> a = np.arange(0, 60, 10).reshape(-1, 1)
>>> a
array([[ 0], [10], [20], [30], [40], [50]])
>>> a.shape
(6, 1)
次に、shapeが(5,)の1次元配列bを作成します:
>>> b = np.arange(0, 5)
>>> b
array([0, 1, 2, 3, 4])
>>> b.shape
(5,)
a と b の和を計算すると、加算表が得られます。これは a と b のすべての要素の組み合わせの和を計算したもので、shape は (6,5) の配列になります:
>>> c = a + b
>>> c
array([[ 0, 1, 2, 3, 4],
[10, 11, 12, 13, 14],
[20, 21, 22, 23, 24],
[30, 31, 32, 33, 34],
[40, 41, 42, 43, 44],
[50, 51, 52, 53, 54]])
>>> c.shape
(6, 5)
a と b の shape の長さ(つまり ndim 属性)が異なるため、ルール 1 に基づいて b の shape を a に合わせる必要があります。したがって、b の shape の前に 1 を追加して (1,5) に補完します。これは次の計算を行ったことに相当します:
>>> b.shape=1,5
>>> b
array([[0, 1, 2, 3, 4]])
このようにして、加算演算の 2 つの入力配列の shape はそれぞれ (6,1) と (1,5) になります。ルール 2 に基づいて、出力配列の各軸の長さは入力配列の各軸の長さの最大値であるため、出力配列の shape は (6,5) になります。
b の第 0 軸の長さは 1 であり、a の第 0 軸の長さは 6 であるため、これらを加算するためには、b の第 0 軸の長さを 6 に拡張する必要があります。これは次のように行われます:
>>> b = b.repeat(6,axis=0)
>>> b
array([[0, 1, 2, 3, 4],
[0, 1, 2, 3, 4],
[0, 1, 2, 3, 4],
[0, 1, 2, 3, 4],
[0, 1, 2, 3, 4],
[0, 1, 2, 3, 4]])
a の第 1 軸の長さは 1 であり、b の第 1 軸の長さは 5 であるため、これらを加算するためには、a の第 1 軸の長さを 5 に拡張する必要があります。これは次のように行われます:
>>> a = a.repeat(5, axis=1)
>>> a
array([[ 0, 0, 0, 0, 0],
[10, 10, 10, 10, 10],
[20, 20, 20, 20, 20],
[30, 30, 30, 30, 30],
[40, 40, 40, 40, 40],
[50, 50, 50, 50, 50]])
上記の処理を経て、a と b は対応する要素を加算することができるようになります。
もちろん、numpy が a+b 演算を実行する際、内部では実際に長さ 1 の軸を repeat 関数で拡張することはありません。このようにすると、メモリを無駄に消費することになります。
このようなブロードキャスト計算は非常に一般的であるため、numpy は上記のような a、b 配列を迅速に生成するための方法を提供しています:ogrid オブジェクト:
>>> x,y = np.ogrid[0:5,0:5]
>>> x
array([[0],
[1],
[2],
[3],
[4]])
>>> y
array([[0, 1, 2, 3, 4]])
ogrid は非常に興味深いオブジェクトで、まるで多次元配列のように、スライスのタプルをインデックスとして使用してアクセスし、ブロードキャスト計算に使用できる配列のセットを返します。スライスのインデックスには 2 つの形式があります:
1. 開始値:終了値:ステップ、np.arange (開始値,終了値,ステップ) に似ています。
2. 開始値:終了値:長さ j、3 番目の引数が虚数の場合、返される配列の長さを示します。np.linspace (開始値,終了値,長さ) に似ています:
>>> x, y = np.ogrid[0:1:4j, 0:1:3j]
>>> x
array([[ 0. ],
[ 0.33333333],
[ 0.66666667],
[ 1. ]])
>>> y
array([[ 0. , 0.5, 1. ]])
ogrid の返り値を使用することで、x、y のグリッド面上の各点の値や、x、y、z のグリッド体上の各点の値を簡単に計算できます。以下は、三次元曲面 x * exp (x2 - y2) を描画するプログラムです:
import numpy as np
from enthought.mayavi import mlab
x, y = np.ogrid[-2:2:20j, -2:2:20j]
z = x * np.exp( - x**2 - y**2)
pl = mlab.surf(x, y, z, warp_scale="auto")
mlab.axes(xlabel='x', ylabel='y', zlabel='z')
mlab.outline(pl)
このプログラムは mayavi の mlab ライブラリを使用して、以下のような 3D 曲面を迅速に描画します。
8. ufunc のメソッド#
reduce メソッドは Python の reduce 関数に似ており、axis 軸に沿って array に対して操作を行います。これは演算子を axis 軸に沿ったすべてのサブ配列または要素の間に挿入することに相当します。
例えば:
>>> np.add.reduce([1,2,3]) # 1 + 2 + 3
6
>>> np.add.reduce([[1,2,3],[4,5,6]], axis=1) # 1,4 + 2,5 + 3,6
array([ 6, 15])
accumulate メソッドは reduce メソッドに似ていますが、返される配列は入力配列と同じ shape を持ち、すべての中間計算結果を保存します:
>>> np.add.accumulate([1,2,3])
array([1, 3, 6])
>>> np.add.accumulate([[1,2,3],[4,5,6]], axis=1)
array([[ 1, 3, 6],
[ 4, 9, 15]])
outer メソッド、.outer (a,b) メソッド:
>>> np.multiply.outer([1,2,3,4,5],[2,3,4])
array([[ 2, 3, 4],
[ 4, 6, 8],
[ 6, 9, 12],
[ 8, 12, 16],
[10, 15, 20]])
outer メソッドで計算された結果は以下のような乗法表です:
# 2, 3, 4
# 1
# 2
# 3
# 4
# 5
もしこれら 2 つの配列を等しいプログラムで一歩ずつ計算した場合、乗法表は最終的にブロードキャスト方式で計算されることになります。
9. 行列運算#
行列の積は dot 関数を使用して計算できます。2 次元配列に対しては行列の積を計算し、1 次元配列に対してはその点積を計算します。1 次元配列を列ベクトルまたは行ベクトルとして行列演算を行う必要がある場合は、最初に reshape 関数を使用して 1 次元配列を 2 次元配列に変換することをお勧めします:
>>> a = array([1, 2, 3])
>>> a.reshape((-1,1))
array([[1],
[2],
[3]])
>>> a.reshape((1,-1))
array([[1, 2, 3]])
dot : 2 つの 1 次元配列に対して、これら 2 つの配列の対応するインデックスの要素の積の和(数学的には内積と呼ばれます)を計算します。実際には点乗算です!
10. ファイルストレージ#
NumPy は配列の内容を保存および取得するためのさまざまなファイル操作関数を提供しています。ファイルの入出力形式は 2 つのカテゴリに分かれます:バイナリとテキスト。バイナリ形式のファイルは、NumPy 専用のフォーマットされたバイナリタイプと無形式タイプに分かれます。
配列のメソッド関数 tofile を使用すると、配列内のデータをバイナリ形式でファイルに簡単に書き込むことができます。tofile で出力されるデータにはフォーマットがないため、numpy.fromfile で読み戻すときは自分でデータをフォーマットする必要があります:
>>> a = np.arange(0,12)
>>> a.shape = 3,4
>>> a
array([[ 0, 1, 2, 3],
[ 4, 5, 6, 7],
[ 8, 9, 10, 11]])
>>> a.tofile("a.bin")
>>> b = np.fromfile("a.bin", dtype=np.float) # float型としてデータを読み込む
>>> b # 読み込まれたデータは誤っています
array([ 2.12199579e-314, 6.36598737e-314, 1.06099790e-313,
1.48539705e-313, 1.90979621e-313, 2.33419537e-313])
>>> a.dtype # aのdtypeを確認
dtype('int32')
>>> b = np.fromfile("a.bin", dtype=np.int32) # int32型としてデータを読み込む
>>> b # データは1次元です
array([ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11])
>>> b.shape = 3, 4 # aのshapeに従ってbのshapeを変更
>>> b # これで正しくなりました
array([[ 0, 1, 2, 3],
[ 4, 5, 6, 7],
[ 8, 9, 10, 11]])
上記の例から、読み込むときに正しい dtype と shape を設定する必要があることがわかります。また、tofile 関数は配列の配置順序が C 言語形式であっても Fortran 言語形式であっても、統一して C 言語形式で出力します。
さらに、fromfile と tofile 関数を呼び出す際に sep キーワードパラメータを指定すると、配列はテキスト形式で入出力されます。
numpy.load と numpy.save 関数は NumPy 専用のバイナリタイプでデータを保存します。これらの 2 つの関数は要素タイプや shape などの情報を自動的に処理するため、これらを使用して配列を読み書きするのは非常に便利ですが、numpy.save で出力されたファイルは他の言語で書かれたプログラムで読み込むのが難しいです:
>>> np.save("a.npy", a)
>>> c = np.load( "a.npy" )
>>> c
array([[ 0, 1, 2, 3],
[ 4, 5, 6, 7],
[ 8, 9, 10, 11]])
複数の配列を 1 つのファイルに保存したい場合は、numpy.savez 関数を使用できます。savez 関数の最初の引数はファイル名で、その後の引数は保存する配列です。また、キーワード引数を使用して配列に名前を付けることもできます。非キーワード引数で渡された配列は自動的に arr_0、arr_1、... という名前が付けられます。savez 関数は圧縮ファイル(拡張子は npz)を出力し、その中の各ファイルは save 関数で保存された npy ファイルであり、ファイル名は配列名に対応します。load 関数は npz ファイルを自動的に認識し、辞書のようなオブジェクトを返します。配列名をキーワードとして使用して配列の内容を取得できます:
>>> a = np.array([[1,2,3],[4,5,6]])
>>> b = np.arange(0, 1.0, 0.1)
>>> c = np.sin(b)
>>> np.savez("result.npz", a, b, sin_array = c)
>>> r = np.load("result.npz")
>>> r["arr_0"] # 配列a
array([[1, 2, 3],
[4, 5, 6]])
>>> r["arr_1"] # 配列b
array([ 0. , 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9])
>>> r["sin_array"] # 配列c
array([ 0. , 0.09983342, 0.19866933, 0.29552021, 0.38941834,
0.47942554, 0.56464247, 0.64421769, 0.71735609, 0.78332691])
もし result.npz ファイルを解凍ソフトで開くと、arr_0.npy、arr_1.npy、sin_array.npy という 3 つのファイルがあり、それぞれに配列 a、b、c の内容が保存されています。
numpy.savetxt と numpy.loadtxt を使用すると、1 次元および 2 次元の配列を読み書きできます:
>>> a = np.arange(0,12,0.5).reshape(4,-1)
>>> np.savetxt("a.txt", a) # デフォルトでは'%.18e'形式でデータを保存し、空白で区切ります
>>> np.loadtxt("a.txt")
array([[ 0. , 0.5, 1. , 1.5, 2. , 2.5],
[ 3. , 3.5, 4. , 4.5, 5. , 5.5],
[ 6. , 6.5, 7. , 7.5, 8. , 8.5],
[ 9. , 9.5, 10. , 10.5, 11. , 11.5]])
>>> np.savetxt("a.txt", a, fmt="%d", delimiter=",") # 整数として保存し、カンマで区切ります
>>> np.loadtxt("a.txt",delimiter=",") # 読み込むときもカンマで区切る必要があります
array([[ 0., 0., 1., 1., 2., 2.],
[ 3., 3., 4., 4., 5., 5.],
[ 6., 6., 7., 7., 8., 8.],
[ 9., 9., 10., 10., 11., 11.]])
この記事は参考にしました: http://old.sebug.net/paper/books/scipydoc/numpy_intro.html