Numbers¶
In [1]:
%matplotlib inline
In [2]:
import matplotlib.pyplot as plt
In [3]:
import numpy as np
The ndarray: Vectors, matrices and tenosrs¶
dtype, shape, strides
Vector¶
In [4]:
x = np.array([1,2,3])
x
Out[4]:
array([1, 2, 3])
In [5]:
type(x)
Out[5]:
numpy.ndarray
In [6]:
x.dtype
Out[6]:
dtype('int64')
In [7]:
x.shape
Out[7]:
(3,)
In [8]:
x.strides
Out[8]:
(8,)
Matrix¶
In [9]:
x = np.array([[1,2,3], [4,5,6]], dtype=np.int32)
x
Out[9]:
array([[1, 2, 3],
[4, 5, 6]], dtype=int32)
In [10]:
x.dtype
Out[10]:
dtype('int32')
In [11]:
x.shape
Out[11]:
(2, 3)
In [12]:
x.strides
Out[12]:
(12, 4)
Tensor¶
In [13]:
x = np.arange(24).reshape((2,3,4))
In [14]:
x
Out[14]:
array([[[ 0, 1, 2, 3],
[ 4, 5, 6, 7],
[ 8, 9, 10, 11]],
[[12, 13, 14, 15],
[16, 17, 18, 19],
[20, 21, 22, 23]]])
Creating ndarrays¶
From a file¶
In [15]:
%%file numbers.txt
a,b,c # can also skip headers
1,2,3
4,5,6
Overwriting numbers.txt
In [16]:
np.loadtxt('numbers.txt', dtype='int', delimiter=',',
skiprows=1, comments='#')
Out[16]:
array([[1, 2, 3],
[4, 5, 6]])
From Python lists or tuples¶
In [17]:
np.array([
[1,2,3],
[4,5,6]
])
Out[17]:
array([[1, 2, 3],
[4, 5, 6]])
From ranges¶
arange, linspace, logspace
In [18]:
np.arange(1, 7).reshape((2,3))
Out[18]:
array([[1, 2, 3],
[4, 5, 6]])
In [19]:
np.linspace(1, 10, 4)
Out[19]:
array([ 1., 4., 7., 10.])
In [20]:
np.logspace(0, 4, 5, dtype='int')
Out[20]:
array([ 1, 10, 100, 1000, 10000])
From a function¶
fromfunciton
In [21]:
np.fromfunction(lambda i, j: i*3 + j + 1, (2,3))
Out[21]:
array([[1., 2., 3.],
[4., 5., 6.]])
In [22]:
np.fromfunction(lambda i, j: (i-2)**2 + (j-2)**2, (5,5), dtype='int')
Out[22]:
array([[8, 5, 4, 5, 8],
[5, 2, 1, 2, 5],
[4, 1, 0, 1, 4],
[5, 2, 1, 2, 5],
[8, 5, 4, 5, 8]])
How to visualize fromfunction¶
In [23]:
j = np.repeat([np.arange(5)], 5, axis=0)
i = j.T
In [24]:
i
Out[24]:
array([[0, 0, 0, 0, 0],
[1, 1, 1, 1, 1],
[2, 2, 2, 2, 2],
[3, 3, 3, 3, 3],
[4, 4, 4, 4, 4]])
In [25]:
j
Out[25]:
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]])
In [26]:
(i-2)**2 + (j-2)**2
Out[26]:
array([[8, 5, 4, 5, 8],
[5, 2, 1, 2, 5],
[4, 1, 0, 1, 4],
[5, 2, 1, 2, 5],
[8, 5, 4, 5, 8]])
Using element-wise functions in fromfunction¶
In [27]:
np.fromfunction(lambda i, j: np.where(i==j,0, -1), (5,5))
Out[27]:
array([[ 0, -1, -1, -1, -1],
[-1, 0, -1, -1, -1],
[-1, -1, 0, -1, -1],
[-1, -1, -1, 0, -1],
[-1, -1, -1, -1, 0]])
In [28]:
np.fromfunction(lambda i, j: np.where(i<j, 1, np.where(i==j,0, -1)), (5,5))
Out[28]:
array([[ 0, 1, 1, 1, 1],
[-1, 0, 1, 1, 1],
[-1, -1, 0, 1, 1],
[-1, -1, -1, 0, 1],
[-1, -1, -1, -1, 0]])
In [29]:
np.fromfunction(lambda i, j: np.minimum(i,j), (5,5), dtype='int')
Out[29]:
array([[0, 0, 0, 0, 0],
[0, 1, 1, 1, 1],
[0, 1, 2, 2, 2],
[0, 1, 2, 3, 3],
[0, 1, 2, 3, 4]])
In [30]:
np.fromfunction(lambda i, j: np.maximum(i,j), (5,5), dtype='int')
Out[30]:
array([[0, 1, 2, 3, 4],
[1, 1, 2, 3, 4],
[2, 2, 2, 3, 4],
[3, 3, 3, 3, 4],
[4, 4, 4, 4, 4]])
From special constructors¶
zeros, ones, eye, diag
In [32]:
np.zeros((2,3))
Out[32]:
array([[0., 0., 0.],
[0., 0., 0.]])
In [33]:
np.ones((2,3))
Out[33]:
array([[1., 1., 1.],
[1., 1., 1.]])
In [34]:
np.eye(3)
Out[34]:
array([[1., 0., 0.],
[0., 1., 0.],
[0., 0., 1.]])
In [35]:
np.eye(3, 4)
Out[35]:
array([[1., 0., 0., 0.],
[0., 1., 0., 0.],
[0., 0., 1., 0.]])
In [36]:
np.eye(4, k=-1)
Out[36]:
array([[0., 0., 0., 0.],
[1., 0., 0., 0.],
[0., 1., 0., 0.],
[0., 0., 1., 0.]])
In [37]:
np.diag([1,2,3,4])
Out[37]:
array([[1, 0, 0, 0],
[0, 2, 0, 0],
[0, 0, 3, 0],
[0, 0, 0, 4]])
In [38]:
np.diag([1,2,3,4], k=1)
Out[38]:
array([[0, 1, 0, 0, 0],
[0, 0, 2, 0, 0],
[0, 0, 0, 3, 0],
[0, 0, 0, 0, 4],
[0, 0, 0, 0, 0]])
From random variables¶
Convenience functions¶
rand, randn
In [39]:
np.random.rand(2,3)
Out[39]:
array([[0.77299768, 0.34178366, 0.23880732],
[0.10302561, 0.5573475 , 0.8210573 ]])
In [40]:
np.random.randn(2,3)
Out[40]:
array([[ 1.71494961, -1.16407052, 0.80896854],
[ 0.8285348 , -1.29637862, 0.91382269]])
Distributions¶
uniform, normal, randint, poisson, multinomial, multivariate_ normal
In [41]:
np.random.uniform(0, 1, (2,3))
Out[41]:
array([[0.4804166 , 0.44907082, 0.04525296],
[0.25025637, 0.06274173, 0.56979373]])
In [42]:
np.random.normal(0, 1, (2,3))
Out[42]:
array([[-0.50581734, -1.75846746, 0.09303513],
[ 0.72487862, 0.80741124, -0.04779211]])
In [43]:
np.random.randint(0, 10, (4,5))
Out[43]:
array([[9, 6, 3, 4, 0],
[6, 9, 7, 0, 9],
[8, 8, 8, 7, 3],
[1, 3, 3, 5, 1]])
In [44]:
np.random.poisson(10, (4,5))
Out[44]:
array([[ 9, 9, 11, 11, 9],
[15, 6, 7, 5, 11],
[ 6, 12, 10, 13, 9],
[ 6, 11, 10, 14, 10]])
In [45]:
np.random.multinomial(n=5, pvals=np.ones(5)/5, size=8)
Out[45]:
array([[1, 0, 2, 0, 2],
[1, 1, 1, 0, 2],
[1, 0, 1, 2, 1],
[1, 1, 1, 2, 0],
[1, 0, 2, 1, 1],
[1, 1, 2, 1, 0],
[1, 0, 1, 2, 1],
[2, 1, 1, 1, 0]])
In [46]:
np.random.multivariate_normal(mean=[10,20,30], cov=np.eye(3), size=4)
Out[46]:
array([[10.17505278, 19.43150301, 28.615404 ],
[ 8.07667886, 19.93391258, 28.61208613],
[ 9.64838835, 21.69728818, 28.87222115],
[10.51988462, 18.76470562, 31.25916781]])
Indexing¶
In [47]:
x = np.arange(20).reshape((4,5))
x
Out[47]:
array([[ 0, 1, 2, 3, 4],
[ 5, 6, 7, 8, 9],
[10, 11, 12, 13, 14],
[15, 16, 17, 18, 19]])
Using slices¶
In [50]:
x[1,:]
Out[50]:
array([5, 6, 7, 8, 9])
In [51]:
x[:,1]
Out[51]:
array([ 1, 6, 11, 16])
In [52]:
x[1:3,1:3]
Out[52]:
array([[ 6, 7],
[11, 12]])
Extrcting blocks with arbitrary row and column lists (fancy indexing)¶
np.ix_
In [54]:
x[:, [0,3]]
Out[54]:
array([[ 0, 3],
[ 5, 8],
[10, 13],
[15, 18]])
Warning: Fancy indexing can only be used for 1 dimension at a time.
In the example below, numpy treats the arguments as paired
coordinates, and returns the values at (0,0) and (2,3).
In [55]:
x[[0,2],[0,3]]
Out[55]:
array([ 0, 13])
Use the helper np.ix_ to extract arbitrary blocks.
In [56]:
x[np.ix_([0,2], [0,3])]
Out[56]:
array([[ 0, 3],
[10, 13]])
A slice is a view, not a copy¶
In [57]:
x
Out[57]:
array([[ 0, 1, 2, 3, 4],
[ 5, 6, 7, 8, 9],
[10, 11, 12, 13, 14],
[15, 16, 17, 18, 19]])
In [58]:
y = x[1:-1, 1:-1]
y
Out[58]:
array([[ 6, 7, 8],
[11, 12, 13]])
In [59]:
y *= 10
In [60]:
y
Out[60]:
array([[ 60, 70, 80],
[110, 120, 130]])
In [61]:
x
Out[61]:
array([[ 0, 1, 2, 3, 4],
[ 5, 60, 70, 80, 9],
[ 10, 110, 120, 130, 14],
[ 15, 16, 17, 18, 19]])
Use the copy method to convert a view to a copy
In [62]:
z = x[1:-1, 1:-1].copy()
In [63]:
z
Out[63]:
array([[ 60, 70, 80],
[110, 120, 130]])
In [64]:
z[:] = 0
In [65]:
z
Out[65]:
array([[0, 0, 0],
[0, 0, 0]])
In [66]:
x
Out[66]:
array([[ 0, 1, 2, 3, 4],
[ 5, 60, 70, 80, 9],
[ 10, 110, 120, 130, 14],
[ 15, 16, 17, 18, 19]])
Boolean indexing¶
In [67]:
x[x % 2 == 0]
Out[67]:
array([ 0, 2, 4, 60, 70, 80, 10, 110, 120, 130, 14, 16, 18])
In [68]:
x [x > 3]
Out[68]:
array([ 4, 5, 60, 70, 80, 9, 10, 110, 120, 130, 14, 15, 16,
17, 18, 19])
Functions that return indexes¶
In [69]:
idx = np.nonzero(x)
idx
Out[69]:
(array([0, 0, 0, 0, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3]),
array([1, 2, 3, 4, 0, 1, 2, 3, 4, 0, 1, 2, 3, 4, 0, 1, 2, 3, 4]))
In [70]:
x[idx]
Out[70]:
array([ 1, 2, 3, 4, 5, 60, 70, 80, 9, 10, 110, 120, 130,
14, 15, 16, 17, 18, 19])
In [71]:
idx = np.where(x > 3)
idx
Out[71]:
(array([0, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3]),
array([4, 0, 1, 2, 3, 4, 0, 1, 2, 3, 4, 0, 1, 2, 3, 4]))
In [72]:
x[idx]
Out[72]:
array([ 4, 5, 60, 70, 80, 9, 10, 110, 120, 130, 14, 15, 16,
17, 18, 19])
Margins and the axis argument¶
In [73]:
x
Out[73]:
array([[ 0, 1, 2, 3, 4],
[ 5, 60, 70, 80, 9],
[ 10, 110, 120, 130, 14],
[ 15, 16, 17, 18, 19]])
The 0th axis has 4 items, the 1st axis has 5 items.
In [74]:
x.shape
Out[74]:
(4, 5)
In [75]:
x.mean()
Out[75]:
35.15
Marginalizing out the 0th axis = column summaries¶
In [76]:
x.mean(axis=0)
Out[76]:
array([ 7.5 , 46.75, 52.25, 57.75, 11.5 ])
Marginalizing out the 1st axis = row summaries¶
In [77]:
x.mean(axis=1)
Out[77]:
array([ 2. , 44.8, 76.8, 17. ])
Note marginalizing out the last axis is a common default.
In [78]:
x.mean(axis=-1)
Out[78]:
array([ 2. , 44.8, 76.8, 17. ])
Marginalization works for higher dimensions in the same way¶
In [79]:
x = np.random.random((2,3,4))
x
Out[79]:
array([[[0.4084202 , 0.99781141, 0.75177799, 0.05699364],
[0.67909501, 0.06607691, 0.04156904, 0.24781812],
[0.16150778, 0.79658604, 0.58497698, 0.77213094]],
[[0.52253736, 0.88397538, 0.80045776, 0.94953923],
[0.24744783, 0.21859536, 0.06607039, 0.73061305],
[0.65047809, 0.79473921, 0.45947565, 0.73221718]]])
In [80]:
x.shape
Out[80]:
(2, 3, 4)
In [81]:
x.mean(axis=0).shape
Out[81]:
(3, 4)
In [82]:
x.mean(axis=1).shape
Out[82]:
(2, 4)
In [83]:
x.mean(axis=2).shape
Out[83]:
(2, 3)
In [84]:
x.mean(axis=(0,1)).shape
Out[84]:
(4,)
In [85]:
x.mean(axis=(0,2)).shape
Out[85]:
(3,)
In [86]:
x.mean(axis=(1,2)).shape
Out[86]:
(2,)
Broadcasting¶
Broadcasting is what happens when numpy tries to perform binary
operations on two arrays with different shapes. In general, shapes are
promoted to make the arrays compatible using the following rule
For each axis from highest to lowest
If both dimensions are the same, do nothing
If one of the dimensions is 1 or None and the other is \(k\), promote to \(k\)
Otherwise print error message
In [87]:
x = np.zeros((3,2))
x.shape
Out[87]:
(3, 2)
In [88]:
x
Out[88]:
array([[0., 0.],
[0., 0.],
[0., 0.]])
Shapes are compatible
In [89]:
y = np.ones(2)
y.shape
Out[89]:
(2,)
In [90]:
x + y
Out[90]:
array([[1., 1.],
[1., 1.],
[1., 1.]])
Shapes are compatible
In [91]:
y = np.ones((1,2))
y.shape
Out[91]:
(1, 2)
In [92]:
x + y
Out[92]:
array([[1., 1.],
[1., 1.],
[1., 1.]])
Shapes are incompatible but can be made compaible by adding empty dimension
In [93]:
y = np.ones(3)
y.shape
Out[93]:
(3,)
In [94]:
try:
x + y
except ValueError as e:
print(e)
operands could not be broadcast together with shapes (3,2) (3,)
In [95]:
y[:, None].shape
Out[95]:
(3, 1)
In [96]:
x + y[:, None]
Out[96]:
array([[1., 1.],
[1., 1.],
[1., 1.]])
Shapes are incompatible
In [97]:
y = np.ones((2,2))
y.shape
Out[97]:
(2, 2)
In [98]:
try:
x + y
except ValueError as e:
print(e)
operands could not be broadcast together with shapes (3,2) (2,2)
More examples of broadcasting¶
In [99]:
x1 = np.arange(12)
In [100]:
x1
Out[100]:
array([ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11])
In [101]:
x1 * 10
Out[101]:
array([ 0, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100, 110])
In [102]:
x2 = np.random.randint(0,10,(3,4))
In [103]:
x2
Out[103]:
array([[7, 5, 4, 6],
[8, 3, 8, 8],
[2, 5, 3, 4]])
In [104]:
x2 * 10
Out[104]:
array([[70, 50, 40, 60],
[80, 30, 80, 80],
[20, 50, 30, 40]])
In [105]:
x2.shape
Out[105]:
(3, 4)
Column-wise broadcasting¶
In [106]:
mu = np.mean(x2, axis=0)
mu.shape
Out[106]:
(4,)
In [107]:
x2 - mu
Out[107]:
array([[ 1.33333333, 0.66666667, -1. , 0. ],
[ 2.33333333, -1.33333333, 3. , 2. ],
[-3.66666667, 0.66666667, -2. , -2. ]])
In [108]:
(x2 - mu).mean(axis=0)
Out[108]:
array([-2.96059473e-16, 2.96059473e-16, 0.00000000e+00, 0.00000000e+00])
Row wise broadcasting¶
In [109]:
mu = np.mean(x2, axis=1)
mu.shape
Out[109]:
(3,)
In [110]:
try:
x2 - mu
except ValueError as e:
print(e)
operands could not be broadcast together with shapes (3,4) (3,)
We can add a “dummy” axis using None or np.newaxis¶
In [111]:
mu[:, None].shape
Out[111]:
(3, 1)
In [112]:
x2 - mu[:, None]
Out[112]:
array([[ 1.5 , -0.5 , -1.5 , 0.5 ],
[ 1.25, -3.75, 1.25, 1.25],
[-1.5 , 1.5 , -0.5 , 0.5 ]])
In [113]:
x2 - mu[:, np.newaxis]
Out[113]:
array([[ 1.5 , -0.5 , -1.5 , 0.5 ],
[ 1.25, -3.75, 1.25, 1.25],
[-1.5 , 1.5 , -0.5 , 0.5 ]])
In [114]:
np.mean(x2 - mu[:, None], axis=1)
Out[114]:
array([0., 0., 0.])
Reshaping works too¶
In [115]:
x2 - mu.reshape((-1,1))
Out[115]:
array([[ 1.5 , -0.5 , -1.5 , 0.5 ],
[ 1.25, -3.75, 1.25, 1.25],
[-1.5 , 1.5 , -0.5 , 0.5 ]])
Exercise in broadcasting¶
Creating a 12 by 12 multiplication table
In [116]:
x = np.arange(1, 13)
x[:,None] * x[None,:]
Out[116]:
array([[ 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12],
[ 2, 4, 6, 8, 10, 12, 14, 16, 18, 20, 22, 24],
[ 3, 6, 9, 12, 15, 18, 21, 24, 27, 30, 33, 36],
[ 4, 8, 12, 16, 20, 24, 28, 32, 36, 40, 44, 48],
[ 5, 10, 15, 20, 25, 30, 35, 40, 45, 50, 55, 60],
[ 6, 12, 18, 24, 30, 36, 42, 48, 54, 60, 66, 72],
[ 7, 14, 21, 28, 35, 42, 49, 56, 63, 70, 77, 84],
[ 8, 16, 24, 32, 40, 48, 56, 64, 72, 80, 88, 96],
[ 9, 18, 27, 36, 45, 54, 63, 72, 81, 90, 99, 108],
[ 10, 20, 30, 40, 50, 60, 70, 80, 90, 100, 110, 120],
[ 11, 22, 33, 44, 55, 66, 77, 88, 99, 110, 121, 132],
[ 12, 24, 36, 48, 60, 72, 84, 96, 108, 120, 132, 144]])
Scaling to have zero mean and unit standard devation for each feature.
In [117]:
x = np.random.normal(10, 5,(3,4))
x
Out[117]:
array([[ 7.11526071, 6.66284037, 13.86245668, 15.28157344],
[17.53107011, 6.92335358, 20.84128098, 6.22506193],
[20.20736019, 9.6887015 , 15.30937277, 12.03098673]])
Scaling column-wise
In [118]:
(x - x.mean(axis=0))/x.std(axis=0)
Out[118]:
array([[-1.38749581, -0.80010784, -0.93389282, 1.09512014],
[ 0.45680586, -0.60983251, 1.38666547, -1.3225013 ],
[ 0.93068995, 1.40994034, -0.45277265, 0.22738116]])
Scaling row-wise
In [119]:
(x - x.mean(axis=1)[:, None])/x.std(axis=1)[:, None]
Out[119]:
array([[-0.93239356, -1.04907463, 0.80773607, 1.17373213],
[ 0.72461462, -0.9280853 , 1.24035099, -1.0368803 ],
[ 1.49421221, -1.17049261, 0.253399 , -0.57711859]])
Combining ndarrays¶
In [120]:
x1 = np.zeros((3,4))
x2 = np.ones((3,5))
x3 = np.eye(4)
In [121]:
x1
Out[121]:
array([[0., 0., 0., 0.],
[0., 0., 0., 0.],
[0., 0., 0., 0.]])
In [122]:
x2
Out[122]:
array([[1., 1., 1., 1., 1.],
[1., 1., 1., 1., 1.],
[1., 1., 1., 1., 1.]])
In [123]:
x3
Out[123]:
array([[1., 0., 0., 0.],
[0., 1., 0., 0.],
[0., 0., 1., 0.],
[0., 0., 0., 1.]])
Binding rows when number of columns is the same¶
In [124]:
np.r_[x1, x3]
Out[124]:
array([[0., 0., 0., 0.],
[0., 0., 0., 0.],
[0., 0., 0., 0.],
[1., 0., 0., 0.],
[0., 1., 0., 0.],
[0., 0., 1., 0.],
[0., 0., 0., 1.]])
Binding columns when number of rows is the same¶
In [125]:
np.c_[x1, x2]
Out[125]:
array([[0., 0., 0., 0., 1., 1., 1., 1., 1.],
[0., 0., 0., 0., 1., 1., 1., 1., 1.],
[0., 0., 0., 0., 1., 1., 1., 1., 1.]])
You can combine more than 2 at a time¶
In [126]:
np.c_[x1, x2, x1]
Out[126]:
array([[0., 0., 0., 0., 1., 1., 1., 1., 1., 0., 0., 0., 0.],
[0., 0., 0., 0., 1., 1., 1., 1., 1., 0., 0., 0., 0.],
[0., 0., 0., 0., 1., 1., 1., 1., 1., 0., 0., 0., 0.]])
Stacking¶
In [127]:
np.vstack([x1, x3])
Out[127]:
array([[0., 0., 0., 0.],
[0., 0., 0., 0.],
[0., 0., 0., 0.],
[1., 0., 0., 0.],
[0., 1., 0., 0.],
[0., 0., 1., 0.],
[0., 0., 0., 1.]])
In [128]:
np.hstack([x1, x2])
Out[128]:
array([[0., 0., 0., 0., 1., 1., 1., 1., 1.],
[0., 0., 0., 0., 1., 1., 1., 1., 1.],
[0., 0., 0., 0., 1., 1., 1., 1., 1.]])
In [129]:
np.dstack([x2, 2*x2, 3*x2])
Out[129]:
array([[[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.]],
[[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.]],
[[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.]]])
Generic stack with axis argument¶
In [130]:
np.stack([x2, 2*x2, 3*x2], axis=0)
Out[130]:
array([[[1., 1., 1., 1., 1.],
[1., 1., 1., 1., 1.],
[1., 1., 1., 1., 1.]],
[[2., 2., 2., 2., 2.],
[2., 2., 2., 2., 2.],
[2., 2., 2., 2., 2.]],
[[3., 3., 3., 3., 3.],
[3., 3., 3., 3., 3.],
[3., 3., 3., 3., 3.]]])
In [131]:
np.stack([x2, 2*x2, 3*x2], axis=1)
Out[131]:
array([[[1., 1., 1., 1., 1.],
[2., 2., 2., 2., 2.],
[3., 3., 3., 3., 3.]],
[[1., 1., 1., 1., 1.],
[2., 2., 2., 2., 2.],
[3., 3., 3., 3., 3.]],
[[1., 1., 1., 1., 1.],
[2., 2., 2., 2., 2.],
[3., 3., 3., 3., 3.]]])
In [132]:
np.stack([x2, 2*x2, 3*x2], axis=2)
Out[132]:
array([[[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.]],
[[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.]],
[[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.],
[1., 2., 3.]]])
Repetition and tiling¶
For a vector¶
In [133]:
x = np.array([1,2,3])
In [134]:
np.repeat(x, 3)
Out[134]:
array([1, 1, 1, 2, 2, 2, 3, 3, 3])
In [135]:
np.tile(x, 3)
Out[135]:
array([1, 2, 3, 1, 2, 3, 1, 2, 3])
For a matrix¶
In [136]:
x = np.arange(6).reshape((2,3))
x
Out[136]:
array([[0, 1, 2],
[3, 4, 5]])
In [137]:
np.repeat(x, 3)
Out[137]:
array([0, 0, 0, 1, 1, 1, 2, 2, 2, 3, 3, 3, 4, 4, 4, 5, 5, 5])
In [138]:
np.repeat(x, 3, axis=0)
Out[138]:
array([[0, 1, 2],
[0, 1, 2],
[0, 1, 2],
[3, 4, 5],
[3, 4, 5],
[3, 4, 5]])
In [139]:
np.repeat(x, 3, axis=1)
Out[139]:
array([[0, 0, 0, 1, 1, 1, 2, 2, 2],
[3, 3, 3, 4, 4, 4, 5, 5, 5]])
In [140]:
np.tile(x, (3,2))
Out[140]:
array([[0, 1, 2, 0, 1, 2],
[3, 4, 5, 3, 4, 5],
[0, 1, 2, 0, 1, 2],
[3, 4, 5, 3, 4, 5],
[0, 1, 2, 0, 1, 2],
[3, 4, 5, 3, 4, 5]])
Splitting ndarrays¶
In [141]:
x = np.arange(32).reshape((4,8))
In [142]:
x
Out[142]:
array([[ 0, 1, 2, 3, 4, 5, 6, 7],
[ 8, 9, 10, 11, 12, 13, 14, 15],
[16, 17, 18, 19, 20, 21, 22, 23],
[24, 25, 26, 27, 28, 29, 30, 31]])
In [143]:
np.split(x, 4)
Out[143]:
[array([[0, 1, 2, 3, 4, 5, 6, 7]]),
array([[ 8, 9, 10, 11, 12, 13, 14, 15]]),
array([[16, 17, 18, 19, 20, 21, 22, 23]]),
array([[24, 25, 26, 27, 28, 29, 30, 31]])]
In [144]:
np.split(x, 4, axis=1)
Out[144]:
[array([[ 0, 1],
[ 8, 9],
[16, 17],
[24, 25]]), array([[ 2, 3],
[10, 11],
[18, 19],
[26, 27]]), array([[ 4, 5],
[12, 13],
[20, 21],
[28, 29]]), array([[ 6, 7],
[14, 15],
[22, 23],
[30, 31]])]
Vectorization¶
Example 1¶
The operators and functions (ufuncs) in Python are vectorized, and will
work element-wise over all entries in an ndarray.
In [145]:
xs = np.zeros(10, dtype='int')
for i in range(10):
xs[i] = i**2
xs
Out[145]:
array([ 0, 1, 4, 9, 16, 25, 36, 49, 64, 81])
In [146]:
xs = np.arange(10)**2
xs
Out[146]:
array([ 0, 1, 4, 9, 16, 25, 36, 49, 64, 81])
Using ufuncs
In [147]:
np.sqrt(xs)
Out[147]:
array([0., 1., 2., 3., 4., 5., 6., 7., 8., 9.])
In [148]:
np.log1p(xs)
Out[148]:
array([0. , 0.69314718, 1.60943791, 2.30258509, 2.83321334,
3.25809654, 3.61091791, 3.91202301, 4.17438727, 4.40671925])
Example 2¶
Scalar product.
In [149]:
n = 10
xs = np.random.rand(n)
ys = np.random.rand(n)
s = 0
for i in range(n):
s += xs[i] * ys[i]
s
Out[149]:
3.239935512385153
In [150]:
np.dot(xs, ys)
Out[150]:
3.239935512385153
In [151]:
xs @ ys
Out[151]:
3.239935512385153
Example 3¶
\begin{align}
y_0 &= \alpha + \beta_1 x_1 + \beta_2 x_2 \\
y_1 &= \alpha + \beta_1 x_1 + \beta_2 x_2 \\
y_2 &= \alpha + \beta_1 x_1 + \beta_2 x_2 \\
\end{align}
In [152]:
m = 3
n = 2
alpha = np.random.rand(1)
betas = np.random.rand(n,1)
xs = np.random.rand(m,n)
In [153]:
alpha
Out[153]:
array([0.01409626])
In [154]:
betas
Out[154]:
array([[0.59629397],
[0.65997692]])
In [155]:
xs
Out[155]:
array([[0.34227866, 0.17384841],
[0.08140117, 0.13100451],
[0.82822641, 0.30764092]])
Using loops¶
In [156]:
ys = np.zeros((m,1))
for i in range(m):
ys[i] = alpha
for j in range(n):
ys[i] += betas[j] * xs[i,j]
ys
Out[156]:
array([[0.3329309 ],
[0.14909524],
[0.71099858]])
Removing inner loop¶
In [157]:
ys = np.zeros((m,1))
for i in range(m):
ys[i] = alpha + xs[i,:].T @ betas
ys
Out[157]:
array([[0.3329309 ],
[0.14909524],
[0.71099858]])
Removing all loops¶
In [158]:
ys = alpha + xs @ betas
ys
Out[158]:
array([[0.3329309 ],
[0.14909524],
[0.71099858]])
Alternative approach¶
The calculaiton with explicit intercepts and coefficients is common in deep learning, where \(\alpha\) is called the bias (\(b\)) and \(\beta\) are called the weights (\(w\)), and each equation is \(y[i] = b + w[i]*x[i]\).
It is common in statisiics to use an augmented matrix in which the first column is all ones, so that all that is needed is a single matrix multiplicaiotn.
In [159]:
X = np.c_[np.ones(m), xs]
X
Out[159]:
array([[1. , 0.34227866, 0.17384841],
[1. , 0.08140117, 0.13100451],
[1. , 0.82822641, 0.30764092]])
In [160]:
alpha
Out[160]:
array([0.01409626])
In [161]:
betas
Out[161]:
array([[0.59629397],
[0.65997692]])
In [162]:
betas_ = np.concatenate([[alpha], betas])
betas_
Out[162]:
array([[0.01409626],
[0.59629397],
[0.65997692]])
In [163]:
ys = X @ betas_
ys
Out[163]:
array([[0.3329309 ],
[0.14909524],
[0.71099858]])
Simulating diffusion¶
In [164]:
w = 100
h = 100
x = np.zeros((w+2,h+2), dtype='float')
x[(w//2-1):(w//2+2), (h//2-1):(h//2+2)] = 1
wts = np.ones(5)/5
for i in range(41):
if i % 10 == 0:
plt.figure()
plt.imshow(x[1:-1, 1:-1], interpolation='nearest')
center = x[1:-1, 1:-1]
left = x[:-2, 1:-1]
right = x[2:, 1:-1]
bottom = x[1:-1, :-2]
top = x[1:-1, 2:]
nbrs = np.dstack([center, left, right, bottom, top])
x = np.sum(wts * nbrs, axis=-1)
In [ ]: