numpy – 行列（ndarray)

2019-10-20 / tau / コメントする

ベクトルと行列の定義

リテラル

ベクトルはnp.array()で引数にリストを指定して定義。

v = np.array([1, 2, 3])
print(v)

# [1 2 3]

v = np.array([1, 2, 3])

print(v)

# [1 2 3]

行列は同じくnp.array()で引数に二次元配列のリストを指定して定義。

m = np.array([[1, 2, 3], [4, 5, 6]])
print(m)

# [[1 2 3]
#  [4 5 6]]

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

print(m)

# [[1 2 3]

# [4 5 6]]

単位行列

numpy.identity(n)でn×nの単位行列を生成。

print(np.identity(4))

# [[1. 0. 0. 0.]
#  [0. 1. 0. 0.]
#  [0. 0. 1. 0.]
#  [0. 0. 0. 1.]]

print(np.identity(4))

# [[1. 0. 0. 0.]

# [0. 1. 0. 0.]

# [0. 0. 1. 0.]

# [0. 0. 0. 1.]]

転置

行列の転置にはtranspose()メソッドを使う。代替として.Tとしてもよい。

m = np.array([[1, 2, 3], [4, 5, 6]])
print(m.transpose())
print(m.T)

# [[1 4]
#  [2 5]
#  [3 6]]
# [[1 4]
#  [2 5]
#  [3 6]]

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

print(m.transpose())

print(m.T)

# [[1 4]

# [2 5]

# [3 6]]

# [[1 4]

# [2 5]

# [3 6]]

一次元配列で定義したベクトルにはtranspose()は効かない。列ベクトルに変換するにはreshape()メソッドを使う（reshape(行数, 列数)）。

v = np.array([1, 2, 3])
print(v.T)
print(v.reshape(3, 1))

# [1 2 3]
# [[1]
#  [2]
#  [3]]

v = np.array([1, 2, 3])

print(v.T)

print(v.reshape(3, 1))

# [1 2 3]

# [[1]

# [2]

# [3]]

演算

定数倍

ベクトル・行列の定数倍は、各要素の定数倍。

a = 2
v = np.array([1, 2, 3])
m = np.array([[1, 2, 3], [4, 5, 6]])

print(a*v)
print(a*m)

# [2 4 6]
# [[ 2  4  6]
#  [ 8 10 12]]

a = 2

v = np.array([1, 2, 3])

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

print(a*v)

print(a*m)

# [2 4 6]

# [[ 2 4 6]

# [ 8 10 12]]

加減

同じ要素数のベクトル、同じ次元の行列同士の下限は要素同士の加減

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

print(a + b)
print(a - b)
print(A + B)
print(A - B)

# [3 6 9]
# [-1 -2 -3]
# [[ 3  6]
#  [ 9 12]]
# [[-1 -2]
#  [-3 -4]]

a = np.array([1, 2, 3])

b = np.array([2, 4, 6])

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

B = np.array([[2, 4], [6, 8]])

print(a + b)

print(a - b)

print(A + B)

print(A - B)

# [3 6 9]

# [-1 -2 -3]

# [[ 3 6]

# [ 9 12]]

# [[-1 -2]

# [-3 -4]]

ベクトルの内積

同じ要素数のベクトルの内積（ドット積）はnp.dot()で計算。

$\begin{equation*} {\bf a} \cdot {\bf b} = \sum_{i=1}^n a_i b_i \end{equation*}$

a = np.array([1, 2, 3])
b = np.array([4, 5, 6])
print(np.dot(a, b))

# 32

a = np.array([1, 2, 3])

b = np.array([4, 5, 6])

print(np.dot(a, b))

# 32

*演算子を使うと、要素ごとの積になる。

print(a*b)

# [ 4 10 18]

print(a*b)

# [ 4 10 18]

行列の積

行列同士の積もnp.dot()で計算。l×m行列とm×n行列の積はl×n行列になる。

$\begin{equation*} \left( \begin{array}{ccc} a_{11} & \cdots & a_{1m} \\ \vdots & a_{ij} & \vdots \\ a_{l1} & \cdots & a_{lm} \\ \end{array} \right) \cdot \left( \begin{array}{ccc} b_{11} & \cdots & b_{1n} \\ \vdots & b_{jk} & \vdots \\ b_{m1} & \cdots & b_{mn} \\ \end{array} \right) = \left( \sum_{j=1}^m a_{ij} b_{jk} \right) \end{equation*}$

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

# [[22 28]
#  [49 64]]

A = np.array([[1, 2, 3],

[4, 5, 6]])

B = np.array([[1, 2],

[3, 4],

[5, 6]])

print(np.dot(A, B))

# [[22 28]

# [49 64]]

次元数が整合しないとエラーになる。

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

# ValueError: shapes (2,3) and (2,2) not aligned: 3 (dim 1) != 2 (dim 0)

A = np.array([[1, 2, 3],

[4, 5, 6]])

B = np.array([[1, 2],

[3, 4]])

print(np.dot(A, B))

# ValueError: shapes (2,3) and (2,2) not aligned: 3 (dim 1) != 2 (dim 0)

行ベクトルと行列の積は、ベクトルを前からかけてok。

$\begin{equation*} (1, 2, 3) \left( \begin{array}{ccc} 1 & 2 & 3 \\ 4 & 5 & 6 \\ 7 & 8 & 9 \end{array} \right) = (30, 36, 42) \end{equation*}$

a = np.array([1, 2, 3])
A = np.array([[1, 2, 3],
              [4, 5, 6],
              [7, 8, 9]])
print(np.dot(a, A))

# [30 36 42]

a = np.array([1, 2, 3])

A = np.array([[1, 2, 3],

[4, 5, 6],

[7, 8, 9]])

print(np.dot(a, A))

# [30 36 42]

行列と列ベクトルの積は、一次元配列のベクトルをreshape()で列ベクトルに変換してから。

$\begin{equation*} \left( \begin{array}{ccc} 1 & 2 & 3 \\ 4 & 5 & 6 \\ 7 & 8 & 9 \end{array} \right) \left( \begin{array}{c} 1 \\ 2 \\ 3 \end{array} \right) = \left( \begin{array}{c} 14 \\ 32 \\ 50 \end{array} \right) \end{equation*}$

A = np.array([[1, 2, 3],
              [4, 5, 6],
              [7, 8, 9]])
a = np.array([1, 2, 3]).reshape(3, 1)
print(np.dot(A, a))

# [[14]
#  [32]
#  [50]]

A = np.array([[1, 2, 3],

[4, 5, 6],

[7, 8, 9]])

a = np.array([1, 2, 3]).reshape(3, 1)

print(np.dot(A, a))

# [[14]

# [32]

# [50]]

なお、np.dot()の代わりに演算子@が使える。ベクトル同士なら内積、少なくともいずれか一つが行列なら行列積。

a = np.array([1, 2])
b = np.array([1, 2, 3]).reshape(3, 1)
A = np.array([[1, 2, 3],
              [4, 5, 6]])
B = np.array([[1, 2],
              [3, 4],
              [5, 6]])
print(a @ a)
print(a @ A)
print(A @ b)
print(A @ B)

# 5
# [ 9 12 15]
# [[14]
#  [32]]
# [[22 28]
#  [49 64]]

a = np.array([1, 2])

b = np.array([1, 2, 3]).reshape(3, 1)

A = np.array([[1, 2, 3],

[4, 5, 6]])

B = np.array([[1, 2],

[3, 4],

[5, 6]])

print(a @ a)

print(a @ A)

print(A @ b)

print(A @ B)

# 5

# [ 9 12 15]

# [[14]

# [32]]

# [[22 28]

# [49 64]]

numpy.linalgパッケージ

行列式

行列式はnumpy.linalgパッケージのdet()関数で得られる。linalgは”linear algebra”の略で、慣例としてLAという名前で代替される。

import numpy.linalg as LA

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

print(LA.det(A))
print(LA.det(B))

# -2.0000000000000004
# 6.000000000000001

import numpy.linalg as LA

A = np.array([[1, 2],

[3, 4]])

B = np.array([[1, 2, 3],

[2, 1, 3],

[3, 1, 2]])

print(LA.det(A))

print(LA.det(B))

# -2.0000000000000004

# 6.000000000000001

逆行列

逆行列はnumpy.linalgパッケージのinv()関数で得られる。

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

print(LA.inv(A))
print(LA.inv(B))

# [[-2.   1. ]
#  [ 1.5 -0.5]]
# [[-0.16666667 -0.16666667  0.5       ]
#  [ 0.83333333 -1.16666667  0.5       ]
#  [-0.16666667  0.83333333 -0.5       ]]

A = np.array([[1, 2],

[3, 4]])

B = np.array([[1, 2, 3],

[2, 1, 3],

[3, 1, 2]])

print(LA.inv(A))

print(LA.inv(B))

# [[-2. 1. ]

# [ 1.5 -0.5]]

# [[-0.16666667 -0.16666667 0.5 ]

# [ 0.83333333 -1.16666667 0.5 ]

# [-0.16666667 0.83333333 -0.5 ]]

固有値・固有ベクトル

正方行列の固有値と固有ベクトルを、eig()関数で得ることができる(行列の固有値・固有ベクトルの例題を用いた)。

A = np.array([[3, 1],
              [2, 4]])
print(LA.eig(A))

# (array([2., 5.]), array([[-0.70710678, -0.4472136 ],
#        [ 0.70710678, -0.89442719]]))

A = np.array([[3, 1],

[2, 4]])

print(LA.eig(A))

# (array([2., 5.]), array([[-0.70710678, -0.4472136 ],

# [ 0.70710678, -0.89442719]]))

結果は固有値が並んだベクトルと固有ベクトルが並んだ配列で、それぞれを取り出して利用する。なお、固有ベクトルはノルムが1となるように正規化されている。

eigenvalues, eigenvectors = LA.eig(A)
print(eigenvalues)
print(eigenvectors)

# [2. 5.]
# [[-0.70710678 -0.4472136 ]
#  [ 0.70710678 -0.89442719]]

eigenvalues, eigenvectors = LA.eig(A)

print(eigenvalues)

print(eigenvectors)

# [2. 5.]

# [[-0.70710678 -0.4472136 ]

# [ 0.70710678 -0.89442719]]

注意が必要なのは固有ベクトルの方で、各固有ベクトルは配列の列ベクトルとして並んでいる。固有ベクトルを取り出す方法は2通り。

固有値に対応するサフィックスで列ベクトルを取り出す。この方法はnumpyの公式ドキュメントにも以下のように書かれている。

v(…, M, M) array: The normalized (unit “length”) eigenvectors, such that the column v[:,i] is the eigenvector corresponding to the eigenvalue w[i].

print(eigenvectors[:,0], eigenvectors[:,1])

# [-0.70710678  0.70710678] [-0.4472136  -0.89442719]

print(eigenvectors[:,0], eigenvectors[:,1])

# [-0.70710678 0.70710678] [-0.4472136 -0.89442719]

固有ベクトルの配列を転置して、行ベクトルの並びにする。

eigenvectors = eigenvectors.T
print(eigenvectors[0], eigenvectors[1])

# [-0.70710678 0.70710678] [-0.4472136 -0.89442719]

eigenvectors = eigenvectors.T

print(eigenvectors[0], eigenvectors[1])

# [-0.70710678 0.70710678] [-0.4472136 -0.89442719]

ハノイの塔～再帰処理

2019-10-05 / tau / コメントする

ハノイの塔とは

ハノイの塔とは、以下のような問題。

3本の棒があり、そのうちの1本に直径が異なる複数の板が直径順に重ねられている
1回の移動で1枚の板を移動するが、移動する板はそれより板がない棒か自身より直径が大きな板がある棒にしか移動できない
全ての板が直径順に、異なる棒に移動し終えたら終了

問題例

たとえば板が3枚の場合は以下のような手順になる。

まず一番下の板を動かすため、

最小の板を他の棒に移動(0→2)
次に2番目の板を別の棒に移動(0→1)
最小の板を2番目の板の上に移動(2→1)
最大の板を空いている棒に移動(0→2)

なお、3つの棒を左から0、1、2としている。

次に2番目に大きい板を動かすため

最小の板を空いている棒に移動(1→0)
2番目の板を最大の板の上に移動(1→2)

最後に1番目の板を2番目の板の上に移動して終了(0→2)

再帰手続き化

これの手続を、次のように一般化する。

まず、n枚の山を3つの軸（org、tmp、dst）のorgからdstに移動する関数を考える。

move(n, org, tmp, dst)

n　動かす山の板の枚数
org　動かす山がある軸
tmp　作業用に使う軸
dst　山を動かす先の軸

この関数は、n枚で受け取った山を軸orgから軸dstへ移動するため、この手続きを次のように再帰的に呼び出す。

Step-1：orgからn-1枚の山をtmpに移動。

move(n-1, org, dst, tmp)

Step-2：orgに残った1枚の板をorgに移動

move(1, org, tmp, dst)

Step-3：tmpに移したn-1枚の山をdstに移動

move(n-1, tmp, org, dst)

ただし山の板の枚数が1の場合はそれ以上再帰呼び出しはせず、板があった軸orgと動かした先の軸dstを記録して復帰。

プログラム化

移動過程を軸番号で表示するだけの処理

コード

移動過程を記録するだけならシンプルで、再帰関数を以下のように定義。一枚の板の移動の場合に移動軸を記録・蓄積し、その結果をコンソールに表示している。

def move(n, org, tmp, dst, solution):
    if n == 1:
        solution.append((org, dst))
        return
    else:
        move(n-1, org, dst, tmp, solution)
        move(1, org, tmp, dst, solution)
        move(n-1, tmp, org, dst, solution)


org, tmp, dst = 0, 1, 2
n = 3
solution = []

move(n, org, tmp, dst, solution)

for movement in solution: print(movement)

def move(n, org, tmp, dst, solution):

if n == 1:

solution.append((org, dst))

return

else:

move(n-1, org, dst, tmp, solution)

move(1, org, tmp, dst, solution)

move(n-1, tmp, org, dst, solution)

org, tmp, dst = 0, 1, 2

n = 3

solution = []

move(n, org, tmp, dst, solution)

for movement in solution: print(movement)

実行結果

板が3枚の場合の結果。左側の軸にある一番上の板を、右側の軸に移動することを意味している。

(0, 2)
(0, 1)
(2, 1)
(0, 2)
(1, 0)
(1, 2)
(0, 2)

(0, 2)

(0, 1)

(2, 1)

(0, 2)

(1, 0)

(1, 2)

(0, 2)

クラス化・画面表示

ハノイの塔を解いていく過程を表示させる。単純なコンソール表示であっても、盤の状態保持やその表示など、コアのアルゴリズムとは別の部分でコードが増える。クラスを作っていくと更にコードが増えるが、全体の見通しはよくなる。

HanoiShaftクラス

盤上の軸のクラス。板を置いたり取り除くほか、任意の位置（高さ）にある板が何か知ることができる。板は直径の小さい順に1, 2, 3,…と整数で表し、軸上の指定位置に板がない場合は0。

class HanoiShaft:
    def __init__(self):
        self.shaft = []

    def disk_at(self, layer):
        return 0 if len(self.shaft) < layer + 1 else self.shaft[layer]

    def put(self, disk):
        self.shaft.append(disk)

    def remove(self):
        return self.shaft.pop()

class HanoiShaft:

def __init__(self):

self.shaft = []

def disk_at(self, layer):

return 0 if len(self.shaft) < layer + 1 else self.shaft[layer]

def put(self, disk):

self.shaft.append(disk)

def remove(self):

return self.shaft.pop()

Hanoiクラス

盤の状態を保持し、問題を解き、その結果を表示するクラス。問題を解く部分と、実態の盤の状態保持や表示は別のものだが、それらをまとめている。

class Hanoi:
    def __init__(self, num_disks):
        self.num_disks = num_disks
        self.solution = []
        self.shafts = [HanoiShaft() for n in range(3)]
        org, tmp, dst = 0, 1, 2
        for disk in range(num_disks, 0, -1):
            self.shafts[org].put(disk)

    def solve(self):
        self.solution = []
        self.move(self.num_disks, 0, 1, 2)

    def move(self, n, org, tmp, dst):
        if n == 1:
            self.solution.append((org, dst))
            return
        else:
            self.move(n-1, org, dst, tmp)
            self.move(1, org, tmp, dst)
            self.move(n-1, tmp, org, dst)

    def display_solution(self):
        self.display_status()
        self.solve()
        for movement in self.solution:
            disk = self.shafts[movement[0]].remove()
            self.shafts[movement[1]].put(disk)
            self.display_status()

    def display_status(self):
        max_width = self.num_disks * 2 - 1
        for layer in range(self.num_disks-1, -1, -1):
            layer_str = ''
            for shaft in self.shafts:
                disk = shaft.disk_at(layer)
                if disk == 0:
                    space_size = (max_width - 1) // 2
                    layer_str += ' '*space_size + '|' + ' '*space_size + ' '
                else:
                    disk_size = disk * 2 - 1
                    space_size = (max_width - disk_size) // 2
                    layer_str += \
                        ' '*space_size + '#'*disk_size + ' '*space_size + ' '
            print(layer_str)
        print('-'*(max_width*3+2))
        print()

class Hanoi:

def __init__(self, num_disks):

self.num_disks = num_disks

self.solution = []

self.shafts = [HanoiShaft() for n in range(3)]

org, tmp, dst = 0, 1, 2

for disk in range(num_disks, 0, -1):

self.shafts[org].put(disk)

def solve(self):

self.solution = []

self.move(self.num_disks, 0, 1, 2)

def move(self, n, org, tmp, dst):

if n == 1:

self.solution.append((org, dst))

return

else:

self.move(n-1, org, dst, tmp)

self.move(1, org, tmp, dst)

self.move(n-1, tmp, org, dst)

def display_solution(self):

self.display_status()

self.solve()

for movement in self.solution:

disk = self.shafts[movement[0]].remove()

self.shafts[movement[1]].put(disk)

self.display_status()

def display_status(self):

max_width = self.num_disks * 2 - 1

for layer in range(self.num_disks-1, -1, -1):

layer_str = ''

for shaft in self.shafts:

disk = shaft.disk_at(layer)

if disk == 0:

space_size = (max_width - 1) // 2

layer_str += ' '*space_size + '|' + ' '*space_size + ' '

else:

disk_size = disk * 2 - 1

space_size = (max_width - disk_size) // 2

layer_str += \

' '*space_size + '#'*disk_size + ' '*space_size + ' '

print(layer_str)

print('-'*(max_width*3+2))

print()

実行部分

hanoi = Hanoi(3)
hanoi.display_solution()

1 2	hanoi = Hanoi(3) hanoi.display_solution()

実行結果

板が3枚の場合の結果。

  #     |     |   
 ###    |     |   
#####   |     |   
-----------------

  |     |     |   
 ###    |     |   
#####   |     #   
-----------------

  |     |     |   
  |     |     |   
#####  ###    #   
-----------------

  |     |     |   
  |     #     |   
#####  ###    |   
-----------------

  |     |     |   
  |     #     |   
  |    ###  ##### 
-----------------

  |     |     |   
  |     |     |   
  #    ###  ##### 
-----------------

  |     |     |   
  |     |    ###  
  #     |   ##### 
-----------------

  |     |     #   
  |     |    ###  
  |     |   ##### 
-----------------

# | |

### | |

##### | |

-----------------

| | |

### | |

##### | #

-----------------

| | |

##### ### #

-----------------

| | |

| # |

##### ### |

-----------------

| | |

| # |

| ### #####

-----------------

| | |

# ### #####

-----------------

| | |

| | ###

# | #####

-----------------

| | #

| | ###

| | #####

-----------------

手続回数

もう一度、以下の再帰手続きを考える。

def move(n, org, tmp, dst, solution):
    if n == 1:
        solution.append((org, dst))
        return
    else:
        move(n-1, org, dst, tmp, solution)
        move(1, org, tmp, dst, solution)
        move(n-1, tmp, org, dst, solution)

def move(n, org, tmp, dst, solution):

if n == 1:

solution.append((org, dst))

return

else:

move(n-1, org, dst, tmp, solution)

move(1, org, tmp, dst, solution)

move(n-1, tmp, org, dst, solution)

n個の山の計算でn-1個の山の関数を2つ呼び出し、それらがまたn-2個の山の関数を2つずつ呼び出す。

計算回数はこれらの和なので、n個の山における関数の呼び出し回数 $N_r(n)$ は以下の通り。

$\begin{equation*} N_r(n) = \sum_{i=1}\n 2^{n-1} = \frac{1 - 2^n}{1 - 2} = 2^n - 1 \end{equation*}$

したがって、階乗計算のように一次的な再帰と異なり、nが1増えるごとに計算量が2倍に増えていくことに注意が必要。

Python3 – リストの初期化

2019-10-03 / tau / コメントする

空のリストの生成

空の場合は添字でアクセスできない。

lst = []
print(len(lst)) # 0
print(lst)      # []
print(lst[0])   # Error: list index out of range

lst = []

print(len(lst)) # 0

print(lst) # []

print(lst[0]) # Error: list index out of range

要素が1個のリストの生成

lst = [2]
print(len(lst)) # 1
print(lst)      # [1]
print(lst[0])   # 2

lst = [2]

print(len(lst)) # 1

print(lst) # [1]

print(lst[0]) # 2

要素が複数個のリストの生成

最も一般的な、要素をリテラルで指定したリスト。

lst = [1, 2]
print(len(lst))       # 2
print(lst)            # [1, 2]
print(lst[0], lst[1]) # 1 2

lst = [1, 2]

print(len(lst)) # 2

print(lst) # [1, 2]

print(lst[0], lst[1]) # 1 2

終わりに”,”をつけても同じ。

lst = [1, 2,]
print(lst) # [1, 2]

1 2	lst = [1, 2,] print(lst) # [1, 2]

要素がない場合はエラー。

lst = [,] # Error: invalid syntax

1	lst = [,] # Error: invalid syntax

多次元のリスト

多次元配列のように定義可能。後ろの方の添字が入れ子の内側のリストの添字。

lst = [[0, 1, 2], [3, 4, 5], [6, 7, 8]]
print(lst[0][1], lst[2][1]) # 1 7

1 2	lst = [[0, 1, 2], [3, 4, 5], [6, 7, 8]] print(lst[0][1], lst[2][1]) # 1 7

異なる要素を持つリスト

異なるタイプの要素を持つことができる。

lst = [1, 'two', [3, 'four']]
print(len(lst))  # 3
print(lst[1])    # two
print(lst[2])    # [3, 'four']
print(lst[1][0]) # t
print(lst[2][1]) # four

lst = [1, 'two', [3, 'four']]

print(len(lst)) # 3

print(lst[1]) # two

print(lst[2]) # [3, 'four']

print(lst[1][0]) # t

print(lst[2][1]) # four

繰り返し指定

*演算子

*演算子で繰り返し数を指定。要素の掛け算ではなく、リストの掛け算のように書く。

lst1 = [0]*3
lst2 = [1, 2]*3
print(lst1) # [0, 0, 0]
print(lst2) # [1, 2, 1, 2, 1, 2]

lst1 = [0]*3

lst2 = [1, 2]*3

print(lst1) # [0, 0, 0]

print(lst2) # [1, 2, 1, 2, 1, 2]

これは+演算子によるリストの結合だと考えるとよい。

lst = [0] + [0] + [0]
print(lst) # [0, 0, 0]

1 2	lst = [0] + [0] + [0] print(lst) # [0, 0, 0]

内包表記

要素が固定の場合。

lst = [0 for n in range(3)]
print(lst) # [0, 0, 0]

1 2	lst = [0 for n in range(3)] print(lst) # [0, 0, 0]

要素を変化させる場合。

lst = [n*2 for n in range(3)]
print(lst) [0, 2, 4]

1 2	lst = [n*2 for n in range(3)] print(lst) [0, 2, 4]

*演算子を使う場合の注意点

*演算子は、同一インスタンスへの参照を複製するだけなので、意図した結果にならない場合がある。

lst = [[0]*2]*3
print(lst)    # [[0, 0], [0, 0], [0, 0]]
lst[0][0] = 1
print(lst)    # [[1, 0], [1, 0], [1, 0]]

lst = [[0]*2]*3

print(lst) # [[0, 0], [0, 0], [0, 0]]

lst[0][0] = 1

print(lst) # [[1, 0], [1, 0], [1, 0]]

[0, 0]が3つ複製されるが、これらは同じインスタンスを指すため、一つの要素の変更が他にも反映されてしまう。

内包表記なら毎回インスタンスが生成されるので、意図した結果になる。

lst = [[0]*2 for n in range(3)]
print(lst)    # [[0, 0], [0, 0], [0, 0]]
lst[0][0] = 1
print(lst)    # [[1, 0], [0, 0], [0, 0]]

lst = [[0]*2 for n in range(3)]

print(lst) # [[0, 0], [0, 0], [0, 0]]

lst[0][0] = 1

print(lst) # [[1, 0], [0, 0], [0, 0]]

単回帰分析～コンビニエンスストア

2019-09-29 / tau / コメントする

概要

Pythonを使った統計計算と図示の練習のため、コンビニエンスストアで単回帰分析をやってみた。

コンビニの店舗数は「商業動態統計年報」の2016年データを使い、説明変数として2015年の国勢調査人口、2018年の国土地理院による都道府県面積、2017年道路統計年報の道路実延長データを使った。

計算コードは以下の通り。

import numpy as np
import matplotlib.pyplot as plt

master_data = np.loadtxt('conv-store-data.csv', delimiter=',',
        skiprows=2, usecols=(1, 2), encoding='utf8')
print(master_data)

y = num_stores = master_data[:,0]
x = population = master_data[:,1]

sxy = np.sum((x - x.mean())*(y - y.mean()))
sxx = np.sum((x - x.mean())**2)
syy = np.sum((y - y.mean())**2)
a = sxy / sxx
b = y.mean() - a*x.mean()
r = sxy / np.sqrt(sxx * syy)
n = len(num_stores)
ESS = n * sxy**2 / sxx
TSS = n * syy
R2 = ESS / TSS

max_population = 15000000
max_stores = 8000
plt.xlim(0, max_population)
plt.ylim(0, max_stores)
plt.xlabel('population')
plt.ylabel('number of stores')

plt.scatter(population, num_stores, label='#stores')

x = np.linspace(0, max_population)
plt.plot(x, a*x + b, color='m')

plt.text(1000000, 7600, 'a = ' + str(a))
plt.text(1000000, 7300, 'b = ' + str(b))
plt.text(1000000, 7000, 'r = '+str(r))
plt.text(1000000, 6700, 'R2 = '+str(R2))

plt.show()

import numpy as np

import matplotlib.pyplot as plt

master_data = np.loadtxt('conv-store-data.csv', delimiter=',',

skiprows=2, usecols=(1, 2), encoding='utf8')

print(master_data)

y = num_stores = master_data[:,0]

x = population = master_data[:,1]

sxy = np.sum((x - x.mean())*(y - y.mean()))

sxx = np.sum((x - x.mean())**2)

syy = np.sum((y - y.mean())**2)

a = sxy / sxx

b = y.mean() - a*x.mean()

r = sxy / np.sqrt(sxx * syy)

n = len(num_stores)

ESS = n * sxy**2 / sxx

TSS = n * syy

R2 = ESS / TSS

max_population = 15000000

max_stores = 8000

plt.xlim(0, max_population)

plt.ylim(0, max_stores)

plt.xlabel('population')

plt.ylabel('number of stores')

plt.scatter(population, num_stores, label='#stores')

x = np.linspace(0, max_population)

plt.plot(x, a*x + b, color='m')

plt.text(1000000, 7600, 'a = ' + str(a))

plt.text(1000000, 7300, 'b = ' + str(b))

plt.text(1000000, 7000, 'r = '+str(r))

plt.text(1000000, 6700, 'R2 = '+str(R2))

plt.show()

人口との関係

コンビニ店舗数と人口の散布図と回帰式を以下に示す。

相関係数が極めて高いのは当然で、やはりコンビニの出店計画には人口ファクターが強く作用していることがわかる。

定数項bが店舗数のオーダーに比べてほぼゼロというのも興味深い。

係数aの値からは5万人弱で1店舗ということになるが、人口10万くらいの都市で2店舗しかないことになり、ちょっと少ないような気がする。もしかすると、一定規模以上の市町村単位や都市単位くらいで層化して出店計画を立てているのかもしれない。

面積との関係

店舗数と面積の関係

次に店舗数と都道府県面積の関係を見てみた。

この結果はかなりはずれで、データを見ても人口が少なく面積が群を抜いて大きい北海道の影響を大きく受けている。

ここで、面積が極端に大きい北海道（83423.83㎢～1位、2906店舗～5位）と面積が小さいが集積度が極端に高い東京都（2193.96㎢～45位、7003店舗～1位）の2つを除いて計算してみる。

これはさらにおかしな結果で、面積が小さいほど店舗数が多いことになる。

考えてみれば、集客数を期待するなら人口が集積している地域が有利だから、人口密度に比例する可能性を考えた方がいいのかもしれない。もし面積が小さい県の方が集積度が高いと想定すると、面積だけを取り出したときに逆の関係になるとも考えられるが、相関係数や決定係数が小さすぎるので考察は難しい。

人口密度

以下は人口密度との関係。

今度はかなりきれいに相関の高さが出ている。

直接的な計算式に入れているかどうかわからないが、GISなどで出店計画を立てるとしたら、人口密度の高いエリアを選んでいくだろうことが想定される。

ただ、店舗数は人口などといった売り上げに直結するデータから導かれるのが普通で、人口密度が高くても人口が少なければ出店インセンティブにはならない。

人口と人口密度の関係

試しに人口を説明変数、人口密度を被説明変数として両者の関係を見てみると、驚くことに「人口が多いほど人口密度が高くなる（あるいはその逆）」という関係になる。

ここから先は人口論や地域論になりそうなので置いておくが、少なくとも日本においては、「狭いところほど人が集まっている傾向がある」ということになりそうである。

もちろんこれは他の国でも一般に当てはまることかもしれないが、朝のラッシュ時に特定の車両に無理やり乗り込んでいる割に離れた車両がすいているとか、1本電車を遅らせたらガラガラだったとか、そのあたりの行動パターンを見ていると、何となく日本に特有のような気がする。

道路延長との関係

最後に道路延長との関係を見てみる。

東京のように稠密な都市は例外とすると、概ね関係はありそうである。ただし相関係数、決定係数は高くない。

コンビニ店舗が道路の利便性に依っていることは推測できるが、やはり人口という売り上げ直結のデータに比べると関係は弱い。

高速道路の延長についても見てみたが、こちらはほとんど関係は見られなかった。

ただ、高速道路の延伸に伴ってコンビニエンスストアの店舗数が伸びているようであり、マクロな延長というよりも物流上のインパクトが大きいことは予想される。

Python3 – 浮動小数点の誤差

2019-09-28 / tau / コメントする

概要

Pythonで浮動小数点誤差に行き当たった。

コンソール上にはErrorが表示され計算結果はnanになるが、配列を使っていたので、問題がどこで生じているかわかるまでにちょっと手間取った。

解析上の考え

中心が(cx, cy)、半径がrの円を考え、x座標を与えてそれに対する円上のy座標を計算しようとしていたとき。

このような円の方程式は、陰関数では以下のようになる。

$\begin{equation*} (x - c_x)^2 + (y - c_y)^2 = r^2 \end{equation*}$

これをy ≥ 0で解くと、

$\begin{equation*} y = c_y + \sqrt{r^2 - (x - c_x)^2} \end{equation*}$

問題の発生

中心の座標が(0.5, 0.5)、半径が0.3となるような円が欲しかったので、xの定義域を0.3～0.8として次のようなコードをPythonで書いた。実際はy < 0の部分も必要で、円周上だけでなく内部の座標をランダムに計算し、それらの相関係数を計算するというものだったが、本質的には変わらない。

import numpy as np

r = 0.3
cx, cy = 0.5, 0.5

x = np.linspace(0.3, 0.8, 11)
y = np.array([cy + np.sqrt(r**2 - (x - cx)**2)])

cor = np.sum((x - x.mean()) * (y - y.mean()) / x.std() / y.std())
print(cor)

import numpy as np

r = 0.3

cx, cy = 0.5, 0.5

x = np.linspace(0.3, 0.8, 11)

y = np.array([cy + np.sqrt(r**2 - (x - cx)**2)])

cor = np.sum((x - x.mean()) * (y - y.mean()) / x.std() / y.std())

print(cor)

その結果、cor=nanとなってしまう。実行時にsqrtの計算で警告が出ている。

C:\Users\tomo\Google �h���C�u\IT_and_Mobile\dev\python\trouble\float.py:7: RuntimeWarning: invalid value encountered in sqrt
  y = np.array([cy + np.sqrt(r**2 - (x - cx)**2)])
nan

C:\Users\tomo\Google �h��C�u\IT_and_Mobile\dev\python\trouble\float.py:7: RuntimeWarning: invalid value encountered in sqrt

y = np.array([cy + np.sqrt(r**2 - (x - cx)**2)])

nan

そこで配列x、yを表示させてみると、yの最後の要素がnanになっている。

print(x)
print(y)

1 2	print(x) print(y)

[0.3  0.35 0.4  0.45 0.5  0.55 0.6  0.65 0.7  0.75 0.8 ]
[[0.7236068  0.75980762 0.78284271 0.79580399 0.8        0.79580399
  0.78284271 0.75980762 0.7236068  0.66583124        nan]]

[0.3 0.35 0.4 0.45 0.5 0.55 0.6 0.65 0.7 0.75 0.8 ]

[[0.7236068 0.75980762 0.78284271 0.79580399 0.8 0.79580399

0.78284271 0.75980762 0.7236068 0.66583124 nan]]

さらに、sqrtの中を表示させてみる。

for xx in x: print(r**2, (xx - cx)**2)

1	for xx in x: print(r2, (xx - cx)2)

0.09 0.04000000000000001
0.09 0.022500000000000006
0.09 0.009999999999999995
0.09 0.0024999999999999988
0.09 0.0
0.09 0.0025000000000000044
0.09 0.010000000000000018
0.09 0.022500000000000006
0.09 0.03999999999999998
0.09 0.0625
0.09 0.09000000000000002

0.09 0.04000000000000001

0.09 0.022500000000000006

0.09 0.009999999999999995

0.09 0.0024999999999999988

0.09 0.0

0.09 0.0025000000000000044

0.09 0.010000000000000018

0.09 0.022500000000000006

0.09 0.03999999999999998

0.09 0.0625

0.09 0.09000000000000002

最後のところで差し引かれる方の値がごくわずかに大きくなっていて、その結果根号の中が負の値となってnanが発生したらしい。

解決策

linspaceの範囲を円周上の上限・下限から少し内側にすることも考えられるが、今回はnumpy.around()を使った。

numpy.around(a, decimals=0)

これは値aを小数点以下decimalsの桁数で丸めてくれるもので、今回は小数点以下10桁目で丸めてみた。

import numpy as np

r = 0.3
cx, cy = 0.5, 0.5

x = np.linspace(0.3, 0.8, 11)
y = np.array([cy + np.sqrt(r**2 - np.around(x - cx, 10)**2)])

cor = np.sum((x - x.mean()) * (y - y.mean()) / x.std() / y.std())
print(cor)

import numpy as np

r = 0.3

cx, cy = 0.5, 0.5

x = np.linspace(0.3, 0.8, 11)

y = np.array([cy + np.sqrt(r**2 - np.around(x - cx, 10)**2)])

cor = np.sum((x - x.mean()) * (y - y.mean()) / x.std() / y.std())

print(cor)

この結果、正常に計算。

-6.625370822144879

1	-6.625370822144879

注意点

解析上値が等しくなる（差し引きゼロになる）ような場合でも、浮動小数点ではごくわずかな誤差のため、エラーとなることがある
numpy.linspace()を使うと上限・下限の間を等分してくれるのでありがたいが、上限・下限自体が他の値と厳密に等しくなければならない場合には、浮動小数点の誤差が生じるため注意が必要

matplotlib

2019-09-23 / tau / コメントする

pyplot～グラフの描画

patches～各種図形

pyplot.imshow()～画像表示

設定

markers
color
cmap～カラーマップ
linestyle

ヒストグラムと関数の重ね合わせ

Python3 – numpy

2019-09-23 / tau / コメントする

全般

numpyのインストール

randomモジュール

基本機能

配列・ndarray

Tips

配列の応用

集計・確率・統計

各種関数

`where()関数`

Numpy.where()は、配列中の要素を検索してそのインデックスを得たり、配列の条件によって指定した配列から要素を選び出すことができる。

numpy – 統計量と組み込み関数の注意点

2019-09-22 / tau / コメントする

概要

Pythonのnumpyで統計量を扱う場合、分散・共分散で注意を要する点がある。Excelの関数でも同様だが、分散・共分散が標本値からそのまま計算した値か、不偏推定量として計算した値か、ということを意識する必要がある。

ここでは、ndarrayから律義に計算した値と、numpyの組み込み関数から計算した値を比較してみる。

標本数・総和・平均

標本として与えた配列の要素数、全要素の和、総和を要素数で除した値で、そのまま母集団の不偏推定量。

import numpy as np

x = np.array([1, 2, 3, 4, 5])
n = len(x)
sum_x = x.sum()
mean_x = sum_x / n
print(n)
print(sum_x)
print(mean_x, x.mean())

# 5
# 15
# 3.0 3.0

import numpy as np

x = np.array([1, 2, 3, 4, 5])

n = len(x)

sum_x = x.sum()

mean_x = sum_x / n

print(n)

print(sum_x)

print(mean_x, x.mean())

# 5

# 15

# 3.0 3.0

分散・標準偏差の食い違い

`var()`、`std()`は標本の分散・標準偏差

分散・共分散については、提供されたメソッドvar()、std()は標本数で除した標本分散・共分散となる。

x = np.array([1, 2, 3, 4, 5])
n = len(x)
var_x = np.sum((x - x.mean())**2) / n
std_x = np.sqrt(var_x)
print(var_x, x.var())
print(std_x, x.std())

# 2.0 2.0
# 1.4142135623730951 1.4142135623730951

x = np.array([1, 2, 3, 4, 5])

n = len(x)

var_x = np.sum((x - x.mean())**2) / n

std_x = np.sqrt(var_x)

print(var_x, x.var())

print(std_x, x.std())

# 2.0 2.0

# 1.4142135623730951 1.4142135623730951

不偏分散やその平方根を求める場合は、var()、std()の引数でddof=1とする。

`cov()`で与えられるのは不偏推定量

numpy.cov()で二つの一次元配列の分散・共分散行列を得ることができる。分散共分散行列は

$\begin{equation*} \left[ \begin{array}{ll} {\rm Var}(X) & {\rm Cov}(XY) \\ {\rm Cov}(XY) & {\rm Var}(Y) \end{array} \right] \end{equation*}$

で定義されるが、以下の計算結果の対角要素は、先の分散の計算結果(2.0)と異なっている。

x = np.array([1, 2, 3, 4, 5])
y = np.array([1, 2, 3, 4, 5])
print(np.cov(x, y))

# [[2.5 2.5]
#  [2.5 2.5]]

x = np.array([1, 2, 3, 4, 5])

y = np.array([1, 2, 3, 4, 5])

print(np.cov(x, y))

# [[2.5 2.5]

# [2.5 2.5]]

これらの値は、分散を標本分散ではなく不偏推定量で計算しているためで、試しに先の計算を不偏分散で計算してみると結果は整合する。

ちなみに標本分散は偏差の二乗和を標本数nで割り、不偏分散はn-1で割る。

x = np.array([1, 2, 3, 4, 5])
n = len(x)
var_x = np.sum((x - x.mean())**2) / (n - 1)
print(var_x)

# 2.5

x = np.array([1, 2, 3, 4, 5])

n = len(x)

var_x = np.sum((x - x.mean())**2) / (n - 1)

print(var_x)

# 2.5

共分散にも注意

上記では分散について確認したが、共分散についても注意が必要。

共分散は一般に以下で定義される。

$\begin{equation*} {\rm Cov}(X, Y) = E((X - \overline{X}) (Y - \overline{Y})) = \frac{1}{n} \sum_{i=1}^{n} (x_i - \overline{x}) (y_i - \overline{y}) \end{equation*}$

これに従って計算した結果。

x = np.array([1, 2, 3, 4, 5])
y = np.array([1, 2, 3, 4, 5])
n = len(x)
cov_xy = np.sum((x - x.mean())*(y - y.mean())) / n
print(cov_xy)

# 2.0

x = np.array([1, 2, 3, 4, 5])

y = np.array([1, 2, 3, 4, 5])

n = len(x)

cov_xy = np.sum((x - x.mean())*(y - y.mean())) / n

print(cov_xy)

# 2.0

この結果は、numpy.cov()で計算した非対角要素2.5と食い違う。

ここでn→n-1として計算しなおすと結果は2.5となり同じ値となる。

共分散用は、たとえば回帰分析の残差変動を計算するときに利用するが、その時には標本の共分散を使うので、注意が必要。

numpyで統計計算をする際、特に共分散については、律義に定義式から計算するのがよさそう。

numpy – 配列の統計計算（二次元配列）

2019-09-22 / tau / コメントする

概要

二次元配列についても、一次元配列と同様な統計関係のメソッド群がある。

二次元についても、同じ機能のメソッドがnumpy、ndarrayのメソッドとして準備されている。

二次元配列の場合は、最小・総和・平均などの計算を全要素／列単位／行単位のいずれで行うかを区別する。具体例として、以下の3×3配列の最小値min()を全要素、列単位、行単位で計算した例を示す。

m = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
print(np.min(m))
print(np.min(m, axis=0))
print(np.min(m, axis=1))

# 1
# [1 2 3]
# [1 4 7]

m = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])

print(np.min(m))

print(np.min(m, axis=0))

print(np.min(m, axis=1))

# 1

# [1 2 3]

# [1 4 7]

引数axisを指定しないと全要素が対象となり、axis=0を指定すると列ごとの計算結果を一次元の配列に、axis=1を指定すると行ごとの計算結果を一次元の配列にして返す。

axis=1の場合、計算は行ごとに行われるが、結果は列ベクトルではなく行ベクトルの配列となる。

また、numpyのメソッドではなくndarrayのメソッドでも同じように使える。

m = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
print(m.min())
print(m.min(axis=0))
print(m.min(axis=1))

# 1
# [1 2 3]
# [1 4 7]

m = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])

print(m.min())

print(m.min(axis=0))

print(m.min(axis=1))

# 1

# [1 2 3]

# [1 4 7]

最小値・最大値

最小値／最大値を返す。

numpy.min(m, axis=0/1)
numpy.max(m, axis=0/1)
m.min(axis=0/1)
m.max(axis=0/1)

結果は同じなので、numpyのメソッドについて実行例を示す。

m = np.array([
    [1, 2, 3],
    [4, 5, 6],
    [7, 8, 9]])
print(np.min(m))         # 1
print(np.min(m, axis=0)) # [1 2 3]
print(np.min(m, axis=1)) # [1 4 7]
print(np.max(m))         # 9
print(np.max(m, axis=0)) # [7 8 9]
print(np.max(m, axis=1)) # [3 6 9]

m = np.array([

[1, 2, 3],

[4, 5, 6],

[7, 8, 9]])

print(np.min(m)) # 1

print(np.min(m, axis=0)) # [1 2 3]

print(np.min(m, axis=1)) # [1 4 7]

print(np.max(m)) # 9

print(np.max(m, axis=0)) # [7 8 9]

print(np.max(m, axis=1)) # [3 6 9]

総和・総積

一次元配列の全要素の和・積を返す。

numpy.sum(m, axis=0/1)
numpy.prod(m, axis=0/1)
m.sum(axis=0/1)
m.prod(axis=0/1)

m = np.array([
    [1, 2, 3],
    [4, 5, 6],
    [7, 8, 9]])
print(np.sum(m))          # 45
print(np.sum(m, axis=0))  # [12 15 18]
print(np.sum(m, axis=1))  # [ 6 15 24]
print(np.prod(m))         # 362880
print(np.prod(m, axis=0)) # [ 28  80 162]
print(np.prod(m, axis=1)) # [  6 120 504]

m = np.array([

[1, 2, 3],

[4, 5, 6],

[7, 8, 9]])

print(np.sum(m)) # 45

print(np.sum(m, axis=0)) # [12 15 18]

print(np.sum(m, axis=1)) # [ 6 15 24]

print(np.prod(m)) # 362880

print(np.prod(m, axis=0)) # [ 28 80 162]

print(np.prod(m, axis=1)) # [ 6 120 504]

平均・分散・標準偏差

一次元配列の要素の平均、分散、標準偏差を返す。分散は標本分散なので、不偏分散が必要な場合はvar()*n/(n-1)とする（ただしn=len(v)）。

numpy.mean(m, axis=0/1)
numpy.var(m, axis=0/1)
numpy.std(m, axis=0/1)
m.mean(axis=0/1)
m.var(axis=0/1)
m.std(axis=0/1)

m = np.array([
    [1, 2, 3],
    [4, 5, 6],
    [7, 8, 9]])
print(np.mean(m))         # 5.0
print(np.mean(m, axis=0)) # [4. 5. 6.]
print(np.mean(m, axis=1)) # [2. 5. 8.]
print(np.var(m))          # 6.666666666666667
print(np.var(m, axis=0))  # [6. 6. 6.]
print(np.var(m, axis=1))  # [0.66666667 0.66666667 0.66666667]
print(np.std(m))          # 2.581988897471611
print(np.std(m, axis=0))  # [2.44948974 2.44948974 2.44948974]
print(np.std(m, axis=1))  # [0.81649658 0.81649658 0.81649658]

m = np.array([

[1, 2, 3],

[4, 5, 6],

[7, 8, 9]])

print(np.mean(m)) # 5.0

print(np.mean(m, axis=0)) # [4. 5. 6.]

print(np.mean(m, axis=1)) # [2. 5. 8.]

print(np.var(m)) # 6.666666666666667

print(np.var(m, axis=0)) # [6. 6. 6.]

print(np.var(m, axis=1)) # [0.66666667 0.66666667 0.66666667]

print(np.std(m)) # 2.581988897471611

print(np.std(m, axis=0)) # [2.44948974 2.44948974 2.44948974]

print(np.std(m, axis=1)) # [0.81649658 0.81649658 0.81649658]

累積和・累積積

一次元配列の要素について先頭から累積して積・和を計算し、それらを要素とする配列を返す。

二次元配列でaxis=0を指定した場合、列ごとに行方向に累積した値が並べられた二次元配列となり、axis=1を指定した場合は行ごとに列方向に累積した値が並べられた二次元配列となる。

numpy.cumsum(m)
numpy.cumprod(m)
m.cumsum()
m.cumprod()

m = np.array([
    [1, 2, 3],
    [4, 5, 6],
    [7, 8, 9]])
print(np.cumsum(m))
print(np.cumsum(m, axis=0))
print(np.cumsum(m, axis=1))
print(np.cumprod(m))
print(np.cumprod(m, axis=0))
print(np.cumprod(m, axis=1))

# [ 1  3  6 10 15 21 28 36 45]
# [[ 1  2  3]
#  [ 5  7  9]
#  [12 15 18]]
# [[ 1  3  6]
#  [ 4  9 15]
#  [ 7 15 24]]

# [     1      2      6     24    120    720   5040  40320 362880]
# [[  1   2   3]
#  [  4  10  18]
#  [ 28  80 162]]
# [[  1   2   6]
#  [  4  20 120]
#  [  7  56 504]]

m = np.array([

[1, 2, 3],

[4, 5, 6],

[7, 8, 9]])

print(np.cumsum(m))

print(np.cumsum(m, axis=0))

print(np.cumsum(m, axis=1))

print(np.cumprod(m))

print(np.cumprod(m, axis=0))

print(np.cumprod(m, axis=1))

# [ 1 3 6 10 15 21 28 36 45]

# [[ 1 2 3]

# [ 5 7 9]

# [12 15 18]]

# [[ 1 3 6]

# [ 4 9 15]

# [ 7 15 24]]

# [ 1 2 6 24 120 720 5040 40320 362880]

# [[ 1 2 3]

# [ 4 10 18]

# [ 28 80 162]]

# [[ 1 2 6]

# [ 4 20 120]

# [ 7 56 504]]

numpy – 配列の統計計算（一次元配列）

2019-09-22 / tau / コメントする

概要

一次元配列について、様々な統計関係の計算をするメソッド群。

同じ機能のメソッドがnumpy、ndarrayのメソッドとして準備されている。例えば一次元配列vについて最小値を求めるメソッドは、numpy.min(v)とv.min()のいずれも同じ結果を返す。

最小値・最大値

一次元配列の要素のうち最小値／最大値を返す。

numpy.min(v)
numpy.max(v)
v.min()
v.max()

結果は同じなので、numpyのメソッドについて実行例を示す。

v = np.array([1, 2, 3, 4])
print(np.min(v))
print(np.max(v))

# 1
# 4

v = np.array([1, 2, 3, 4])

print(np.min(v))

print(np.max(v))

# 1

# 4

総和・総積

一次元配列の全要素の和・積を返す。

numpy.sum(v)
numpy.prod(v)
v.sum()
v.prod()

v = np.array([1, 2, 3, 4])
print(np.sum(v))
print(np.prod(v))

# 10
# 24

v = np.array([1, 2, 3, 4])

print(np.sum(v))

print(np.prod(v))

# 10

# 24

平均・分散・標準偏差

一次元配列の要素の平均、分散、標準偏差を返す。分散は標本分散なので、不偏分散が必要な場合はvar()*n/(n-1)とする（ただしn=len(v)）。

numpy.mean(v)
numpy.var(v)
numpy.std(v)
v.mean()
v.var()
v.std()

v = np.array([1, 2, 3, 4])
print(np.mean(v))
print(np.var(v))
print(np.std(v))

# 2.5
# 1.25
# 1.118033988749895

v = np.array([1, 2, 3, 4])

print(np.mean(v))

print(np.var(v))

print(np.std(v))

# 2.5

# 1.25

# 1.118033988749895

$\begin{eqnarray*} \overline{V} &=& \frac{1 + 2 + 3 + 4}{4} = 2.5 \\ {\rm Var}(V) &=& \frac{(1-2.5)^2+(2-2.5)^2+(3-2.5)^2+(4-2.5)^2}{4} \\ &=& \frac{2.25+0.25+0.25+2.25}{4} = \frac{5}{4} \\ &=& 1.25 \\ \sigma_V &=& \sqrt{1.25} = 1.1180 \ldots \end{eqnarray*}$

累積和・累積積

一次元配列の要素について先頭から累積して積・和を計算し、それらを要素とする配列を返す。

numpy.cumsum(v)
numpy.cumprod(v)
v.cumsum()
v.cumprod()

v = np.array([1, 2, 3, 4])
print(np.cumsum(v))
print(np.cumprod(v))

# [ 1  3  6 10]
# [ 1  2  6 24]

ベクトルと行列の定義

リテラル

単位行列

転置

演算

定数倍

加減

ベクトルの内積

行列の積

numpy.linalgパッケージ

行列式

逆行列

固有値・固有ベクトル

ハノイの塔とは

問題例

再帰手続き化

プログラム化

移動過程を軸番号で表示するだけの処理

コード

実行結果

クラス化・画面表示

HanoiShaftクラス

Hanoiクラス

実行部分

実行結果

手続回数

空のリストの生成

要素が1個のリストの生成

要素が複数個のリストの生成

多次元のリスト

異なる要素を持つリスト

繰り返し指定

*演算子

内包表記

*演算子を使う場合の注意点

概要

人口との関係

面積との関係

店舗数と面積の関係

人口密度

人口と人口密度の関係

道路延長との関係

概要

解析上の考え

問題の発生

解決策

注意点

全般

基本機能

配列・ndarray

集計・確率・統計

各種関数

where()関数

概要

標本数・総和・平均

分散・標準偏差の食い違い

var()、std()は標本の分散・標準偏差

cov()で与えられるのは不偏推定量

共分散にも注意

概要

最小値・最大値

総和・総積

平均・分散・標準偏差

累積和・累積積

概要

最小値・最大値

総和・総積

平均・分散・標準偏差

累積和・累積積

`where()関数`

`var()`、`std()`は標本の分散・標準偏差

`cov()`で与えられるのは不偏推定量