# Keypoints in this Section

- slice is a view but not a copy

- using `axis` to calculate aggregate function on margins
    - using dummy axis `a[:, None]` to expand the vector


- Generating Arrays
    - `np.fromfunction()`
    
        - `lambda i, j: 2*i+j, (2, 3)`
        
        - `lambda i, j: np.where(i<j, 1, np.where(i==j, 0, -1)`
        
    - `np.eye(3, 4, k = -1)`, `np.diag()`
    
    - `np.random.rand()`, `np.random,randn()`

- Expanding Arrays
    - Use `np.c_[]` and `np.r_` to combine arrays
    
    - `np.stack`: see [SciPy.org](https://docs.scipy.org/doc/numpy-1.15.0/reference/generated/numpy.stack.html)

    - `np.repeat()` and `np.tile()`
    
- Dividng Arrays
    - `np.split()`



# Numbers


```python
%matplotlib inline
```


```python
import matplotlib.pyplot as plt
```


```python
import numpy as np
```

## The `ndarray`: Vectors, matrices and tenosrs

dtype, shape, strides

### Vector


```python
x = np.array([1,2,3])
x
```




    array([1, 2, 3])




```python
type(x)
```




    numpy.ndarray




```python
x.dtype
```




    dtype('int64')




```python
x.shape
```




    (3,)




```python
x.strides
```




    (8,)



### Matrix


```python
x = np.array([[1,2,3], [4,5,6]], dtype=np.int32)
x
```




    array([[1, 2, 3],
           [4, 5, 6]], dtype=int32)




```python
x.dtype
```




    dtype('int32')




```python
x.shape
```




    (2, 3)




```python
x.strides
```




    (12, 4)



### Tensor


```python
x = np.arange(24).reshape((2,3,4))
```


```python
x
```




    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 `ndarray`s

### From a file


```python
%%file numbers.txt
a,b,c # can also skip headers
1,2,3
4,5,6
```

    Overwriting numbers.txt



```python
np.loadtxt('numbers.txt', dtype='int', delimiter=',',
           skiprows=1, comments='#')
```




    array([[1, 2, 3],
           [4, 5, 6]])



### From Python lists or tuples


```python
np.array([
    [1,2,3],
    [4,5,6]
])
```




    array([[1, 2, 3],
           [4, 5, 6]])



### From ranges

arange, linspace, logspace


```python
np.arange(1, 7).reshape((2,3))
```




    array([[1, 2, 3],
           [4, 5, 6]])




```python
np.linspace(1, 10, 4)
```




    array([ 1.,  4.,  7., 10.])




```python
np.logspace(0, 4, 5, dtype='int')
```




    array([    1,    10,   100,  1000, 10000])



### From a function

`fromfunciton`


```python
np.fromfunction(lambda i, j: i*3 + j + 1, (2,3))
```




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




```python
np.fromfunction(lambda i, j: (i-2)**2 + (j-2)**2, (5,5), dtype='int')
```




    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` 


```python
j = np.repeat([np.arange(5)], 5, axis=0)
i = j.T
```


```python
i
```




    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]])




```python
j
```




    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]])




```python
(i-2)**2 + (j-2)**2
```




    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`


```python
np.fromfunction(lambda i, j: np.where(i==j,0, -1), (5,5))
```




    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]])




```python
np.fromfunction(lambda i, j: np.where(i<j, 1, np.where(i==j,0, -1)), (5,5))
```




    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]])




```python
np.fromfunction(lambda i, j: np.minimum(i,j), (5,5), dtype='int')
```




    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]])




```python
np.fromfunction(lambda i, j: np.maximum(i,j), (5,5), dtype='int')
```




    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


```python
np.zeros((2,3))
```




    array([[0., 0., 0.],
           [0., 0., 0.]])




```python
np.ones((2,3))
```




    array([[1., 1., 1.],
           [1., 1., 1.]])




```python
np.eye(3)
```




    array([[1., 0., 0.],
           [0., 1., 0.],
           [0., 0., 1.]])




```python
np.eye(3, 4)
```




    array([[1., 0., 0., 0.],
           [0., 1., 0., 0.],
           [0., 0., 1., 0.]])




```python
np.eye(4, k=-1)
```




    array([[0., 0., 0., 0.],
           [1., 0., 0., 0.],
           [0., 1., 0., 0.],
           [0., 0., 1., 0.]])




```python
np.diag([1,2,3,4])
```




    array([[1, 0, 0, 0],
           [0, 2, 0, 0],
           [0, 0, 3, 0],
           [0, 0, 0, 4]])




```python
np.diag([1,2,3,4], k=1)
```




    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


```python
np.random.rand(2,3)
```




    array([[0.00460297, 0.6790471 , 0.66693293],
           [0.46481665, 0.68104561, 0.47717252]])




```python
np.random.randn(2,3)
```




    array([[ 0.32728775, -0.7593871 ,  2.00243124],
           [-0.45911955, -2.29433348, -0.62233783]])



#### Distributions

uniform, normal, randint, poisson, multinomial, multivariate_ normal


```python
np.random.uniform(0, 1, (2,3))
```




    array([[0.35471781, 0.68890474, 0.69249992],
           [0.89704091, 0.84881458, 0.55436389]])




```python
np.random.normal(0, 1, (2,3))
```




    array([[ 0.38938858,  0.02046035,  1.12586858],
           [ 1.55033516,  1.42547904, -1.02902684]])




```python
np.random.randint(0, 10, (4,5))
```




    array([[3, 8, 3, 9, 6],
           [9, 7, 9, 0, 2],
           [5, 3, 0, 4, 4],
           [8, 2, 7, 3, 0]])




```python
np.random.poisson(10, (4,5))
```




    array([[12,  5, 11,  8, 11],
           [ 8, 13, 12,  9, 16],
           [10, 16, 17,  7, 14],
           [15,  5, 10, 10,  6]])




```python
np.random.multinomial(n=5, pvals=np.ones(5)/5, size=8)
```




    array([[2, 1, 0, 1, 1],
           [1, 1, 2, 0, 1],
           [1, 1, 1, 1, 1],
           [1, 2, 0, 1, 1],
           [1, 1, 1, 1, 1],
           [0, 0, 3, 0, 2],
           [0, 1, 2, 1, 1],
           [1, 3, 0, 1, 0]])




```python
np.random.multivariate_normal(mean=[10,20,30], cov=np.eye(3), size=4)
```




    array([[ 8.46959709, 19.70106618, 30.96751184],
           [ 7.750921  , 21.15980554, 30.28468523],
           [ 8.22849672, 20.22017755, 29.69429188],
           [ 9.715338  , 18.56547751, 29.56025953]])



## Indexing 


```python
x = np.arange(20).reshape((4,5))
x
```




    array([[ 0,  1,  2,  3,  4],
           [ 5,  6,  7,  8,  9],
           [10, 11, 12, 13, 14],
           [15, 16, 17, 18, 19]])



### Extracing a scalar


```python
x[1,1]
```




    6



### Extracting a vector


```python
x[1]
```




    array([5, 6, 7, 8, 9])



### Using slices


```python
x[1,:]
```




    array([5, 6, 7, 8, 9])




```python
x[:,1]
```




    array([ 1,  6, 11, 16])




```python
x[1:3,1:3]
```




    array([[ 6,  7],
           [11, 12]])



### Using slices with strides


```python
x[::2,::2]
```




    array([[ 0,  2,  4],
           [10, 12, 14]])



### Extrcting blocks with arbitrary row and column lists (fancy indexing)

`np.ix_`


```python
x[:, [0,3]]
```




    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).


```python
x[[0,2],[0,3]]
```




    array([ 0, 13])



Use the helper `np.ix_` to extract arbitrary blocks.


```python
x[np.ix_([0,2], [0,3])]
```




    array([[ 0,  3],
           [10, 13]])



### A slice is a view, not a copy


```python
x
```




    array([[ 0,  1,  2,  3,  4],
           [ 5,  6,  7,  8,  9],
           [10, 11, 12, 13, 14],
           [15, 16, 17, 18, 19]])




```python
y = x[1:-1, 1:-1]
y
```




    array([[ 6,  7,  8],
           [11, 12, 13]])




```python
y *= 10
```


```python
y
```




    array([[ 60,  70,  80],
           [110, 120, 130]])




```python
x
```




    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


```python
z = x[1:-1, 1:-1].copy()
```


```python
z
```




    array([[ 60,  70,  80],
           [110, 120, 130]])




```python
z[:] = 0
```


```python
z
```




    array([[0, 0, 0],
           [0, 0, 0]])




```python
x
```




    array([[  0,   1,   2,   3,   4],
           [  5,  60,  70,  80,   9],
           [ 10, 110, 120, 130,  14],
           [ 15,  16,  17,  18,  19]])



### Boolean indexing


```python
x[x % 2 == 0]
```




    array([  0,   2,   4,  60,  70,  80,  10, 110, 120, 130,  14,  16,  18])




```python
x [x > 3]
```




    array([  4,   5,  60,  70,  80,   9,  10, 110, 120, 130,  14,  15,  16,
            17,  18,  19])




```python
x > 3
```




    array([[False, False, False, False,  True],
           [ True,  True,  True,  True,  True],
           [ True,  True,  True,  True,  True],
           [ True,  True,  True,  True,  True]])



### Functions that return indexes


```python
idx = np.nonzero(x)
idx
```




    (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]))




```python
x[idx]
```




    array([  1,   2,   3,   4,   5,  60,  70,  80,   9,  10, 110, 120, 130,
            14,  15,  16,  17,  18,  19])




```python
idx = np.where(x > 3)
idx
```




    (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]))




```python
x[idx]
```




    array([  4,   5,  60,  70,  80,   9,  10, 110, 120, 130,  14,  15,  16,
            17,  18,  19])



## Margins and the `axis` argument


```python
x
```




    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.


```python
x.shape
```




    (4, 5)




```python
x.mean()
```




    35.15



### Marginalizing out the 0th axis = column summaries


```python
x.mean(axis=0)
```




    array([ 7.5 , 46.75, 52.25, 57.75, 11.5 ])



### Marginalizing out the 1st axis = row summaries


```python
x.mean(axis=1)
```




    array([ 2. , 44.8, 76.8, 17. ])



Note marginalizing out the last axis is a common default.


```python
x.mean(axis=-1)
```




    array([ 2. , 44.8, 76.8, 17. ])



### Marginalization works for higher dimensions in the same way


```python
x = np.random.random((2,3,4))
x
```




    array([[[0.18658517, 0.64047713, 0.94587991, 0.61495165],
            [0.05194638, 0.25243009, 0.77169445, 0.95496594],
            [0.0551977 , 0.85674614, 0.92755212, 0.86043415]],
    
           [[0.84571303, 0.89371063, 0.23284013, 0.09526368],
            [0.91167377, 0.27580046, 0.93419312, 0.79156371],
            [0.45611741, 0.95214143, 0.42178673, 0.84316195]]])




```python
x.shape
```




    (2, 3, 4)




```python
x.mean(axis=0).shape
```




    (3, 4)




```python
x.mean(axis=1).shape
```




    (2, 4)




```python
x.mean(axis=2).shape
```




    (2, 3)




```python
x.mean(axis=(0,1)).shape
```




    (4,)




```python
x.mean(axis=(0,2)).shape
```




    (3,)




```python
x.mean(axis=(1,2)).shape
```




    (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


```python
x = np.zeros((3,2))
x.shape
```




    (3, 2)




```python
x
```




    array([[0., 0.],
           [0., 0.],
           [0., 0.]])



Shapes are compatible


```python
y = np.ones(2)
y.shape
```




    (2,)




```python
x + y
```




    array([[1., 1.],
           [1., 1.],
           [1., 1.]])



Shapes are compatible


```python
y = np.ones((1,2))
y.shape
```




    (1, 2)




```python
x + y
```




    array([[1., 1.],
           [1., 1.],
           [1., 1.]])



Shapes are incompatible but can be made compaible by adding empty dimension


```python
y = np.ones(3)
y.shape
```




    (3,)




```python
try:
    x + y
except ValueError as e:
    print(e)
```

    operands could not be broadcast together with shapes (3,2) (3,) 



```python
y[:, None].shape
```




    (3, 1)




```python
x + y[:, None]
```




    array([[1., 1.],
           [1., 1.],
           [1., 1.]])



Shapes are incompatible


```python
y = np.ones((2,2))
y.shape
```




    (2, 2)




```python
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


```python
x1 = np.arange(12)
```


```python
x1
```




    array([ 0,  1,  2,  3,  4,  5,  6,  7,  8,  9, 10, 11])




```python
x1 * 10
```




    array([  0,  10,  20,  30,  40,  50,  60,  70,  80,  90, 100, 110])




```python
x2 = np.random.randint(0,10,(3,4))
```


```python
x2
```




    array([[3, 3, 6, 2],
           [5, 1, 6, 2],
           [8, 0, 2, 8]])




```python
x2 * 10
```




    array([[30, 30, 60, 20],
           [50, 10, 60, 20],
           [80,  0, 20, 80]])




```python
x2.shape
```




    (3, 4)



### Column-wise broadcasting


```python
mu = np.mean(x2, axis=0)
mu.shape
```




    (4,)




```python
x2 - mu
```




    array([[-2.33333333,  1.66666667,  1.33333333, -2.        ],
           [-0.33333333, -0.33333333,  1.33333333, -2.        ],
           [ 2.66666667, -1.33333333, -2.66666667,  4.        ]])




```python
(x2 - mu).mean(axis=0)
```




    array([ 2.96059473e-16,  7.40148683e-17, -2.96059473e-16,  0.00000000e+00])



### Row wise broadcasting


```python
mu = np.mean(x2, axis=1)
mu.shape
```




    (3,)




```python
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`


```python
mu[:, None].shape
```




    (3, 1)




```python
x2 - mu[:, None]
```




    array([[-0.5, -0.5,  2.5, -1.5],
           [ 1.5, -2.5,  2.5, -1.5],
           [ 3.5, -4.5, -2.5,  3.5]])




```python
x2 - mu[:, np.newaxis]
```




    array([[-0.5, -0.5,  2.5, -1.5],
           [ 1.5, -2.5,  2.5, -1.5],
           [ 3.5, -4.5, -2.5,  3.5]])




```python
np.mean(x2 - mu[:, None], axis=1)
```




    array([0., 0., 0.])



#### Reshaping works too


```python
x2 - mu.reshape((-1,1))
```




    array([[-0.5, -0.5,  2.5, -1.5],
           [ 1.5, -2.5,  2.5, -1.5],
           [ 3.5, -4.5, -2.5,  3.5]])



#### Exercise in broadcasting

Creating a 12 by 12 multiplication table


```python
x = np.arange(1, 13)
x[:,None] * x[None,:]
```




    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.


```python
x = np.random.normal(10, 5,(3,4))
x
```




    array([[10.68728899, 14.63808698, 12.13188671, 13.48607715],
           [ 6.3118309 , 11.92912759,  4.01779711,  7.28426382],
           [ 9.15695838, 20.47190356, 18.5117833 ,  3.00547023]])



Scaling column-wise


```python
(x - x.mean(axis=0))/x.std(axis=0)
```




    array([[ 1.08584868, -0.2922203 ,  0.09746085,  1.29242016],
           [-1.32758954, -1.05220344, -1.27056349, -0.14898015],
           [ 0.24174086,  1.34442375,  1.17310264, -1.14344001]])



Scaling row-wise


```python
(x - x.mean(axis=1)[:, None])/x.std(axis=1)[:, None]
```




    array([[-1.38564641,  1.28669248, -0.40851352,  0.50746745],
           [-0.37305115,  1.57824065, -1.16993435, -0.03525515],
           [-0.51242021,  1.08501579,  0.80828739, -1.38088297]])



## Combining `ndarray`s


```python
x1 = np.zeros((3,4))
x2 = np.ones((3,5))
x3 = np.eye(4)
```


```python
x1
```




    array([[0., 0., 0., 0.],
           [0., 0., 0., 0.],
           [0., 0., 0., 0.]])




```python
x2
```




    array([[1., 1., 1., 1., 1.],
           [1., 1., 1., 1., 1.],
           [1., 1., 1., 1., 1.]])




```python
x3
```




    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


```python
np.r_[x1, x3]
```




    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


```python
np.c_[x1, x2]
```




    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


```python
np.c_[x1, x2, x1]
```




    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


```python
np.vstack([x1, x3])
```




    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.]])




```python
np.hstack([x1, x2])
```




    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.]])




```python
np.dstack([x2, 2*x2, 3*x2])
```




    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


```python
np.stack([x2, 2*x2, 3*x2], axis=0)
```




    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.]]])




```python
np.stack([x2, 2*x2, 3*x2], axis=1)
```




    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.]]])




```python
np.stack([x2, 2*x2, 3*x2], axis=2)
```




    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


```python
x = np.array([1,2,3])
```


```python
np.repeat(x, 3)
```




    array([1, 1, 1, 2, 2, 2, 3, 3, 3])




```python
np.tile(x, 3)
```




    array([1, 2, 3, 1, 2, 3, 1, 2, 3])



#### For a matrix


```python
x = np.arange(6).reshape((2,3))
x
```




    array([[0, 1, 2],
           [3, 4, 5]])




```python
np.repeat(x, 3)
```




    array([0, 0, 0, 1, 1, 1, 2, 2, 2, 3, 3, 3, 4, 4, 4, 5, 5, 5])




```python
np.repeat(x, 3, axis=0)
```




    array([[0, 1, 2],
           [0, 1, 2],
           [0, 1, 2],
           [3, 4, 5],
           [3, 4, 5],
           [3, 4, 5]])




```python
np.repeat(x, 3, axis=1)
```




    array([[0, 0, 0, 1, 1, 1, 2, 2, 2],
           [3, 3, 3, 4, 4, 4, 5, 5, 5]])




```python
np.tile(x, (3,2))
```




    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 `ndarray`s


```python
x = np.arange(32).reshape((4,8))
```


```python
x
```




    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]])




```python
np.split(x, 4)
```




    [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]])]




```python
np.split(x, 4, axis=1)
```




    [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`.


```python
xs = np.zeros(10, dtype='int')
for i in range(10):
    xs[i] = i**2
xs
```




    array([ 0,  1,  4,  9, 16, 25, 36, 49, 64, 81])




```python
xs = np.arange(10)**2
xs
```




    array([ 0,  1,  4,  9, 16, 25, 36, 49, 64, 81])



Using ufuncs


```python
np.sqrt(xs)
```




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




```python
np.log1p(xs)
```




    array([0.        , 0.69314718, 1.60943791, 2.30258509, 2.83321334,
           3.25809654, 3.61091791, 3.91202301, 4.17438727, 4.40671925])



### Example 2

Scalar product.


```python
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
```




    2.5188405619306895




```python
np.dot(xs, ys)
```




    2.51884056193069




```python
xs @ ys
```




    2.51884056193069



### 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}





```python
m = 3
n = 2

alpha = np.random.rand(1)
betas = np.random.rand(n,1)
xs = np.random.rand(m,n)
```


```python
alpha
```




    array([0.14766101])




```python
betas
```




    array([[0.73667336],
           [0.0346741 ]])




```python
xs
```




    array([[0.03376605, 0.60183824],
           [0.70389785, 0.04904536],
           [0.22335997, 0.74175809]])



### Using loops


```python
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
```




    array([[0.19340376],
           [0.6679044 ],
           [0.33792414]])



### Removing inner loop


```python
ys = np.zeros((m,1))
for i in range(m):
    ys[i] = alpha + xs[i,:].T @ betas
ys
```




    array([[0.19340376],
           [0.6679044 ],
           [0.33792414]])



### Removing all loops


```python
ys = alpha + xs @ betas
ys
```




    array([[0.19340376],
           [0.6679044 ],
           [0.33792414]])



### 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.


```python
X = np.c_[np.ones(m), xs]
X
```




    array([[1.        , 0.03376605, 0.60183824],
           [1.        , 0.70389785, 0.04904536],
           [1.        , 0.22335997, 0.74175809]])




```python
alpha
```




    array([0.14766101])




```python
betas
```




    array([[0.73667336],
           [0.0346741 ]])




```python
betas_ = np.concatenate([[alpha], betas])
betas_
```




    array([[0.14766101],
           [0.73667336],
           [0.0346741 ]])




```python
ys = X @ betas_
ys
```




    array([[0.19340376],
           [0.6679044 ],
           [0.33792414]])



### Simulating diffusion


```python
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)
```

9fefd5dc3d2a8e9257ab5ab18691751635232e54
ipynb
Jupyter Notebook
notebook/S03_Numpy_Annotated_more.ipynb
ywang512/sta663
722c73a9b549b3bbbc1ec2ddd8d7ba4750023b92
BSD-3-Clause
2020-02-18T10:06:51.000Z
2020-02-18T10:06:51.000Z
notebook/S03_Numpy_Annotated_more.ipynb
ywang512/sta663
722c73a9b549b3bbbc1ec2ddd8d7ba4750023b92
BSD-3-Clause
notebook/S03_Numpy_Annotated_more.ipynb
ywang512/sta663
722c73a9b549b3bbbc1ec2ddd8d7ba4750023b92
BSD-3-Clause
2021-03-15T08:38:14.000Z
2021-03-15T08:38:14.000Z
Qwen/Qwen-72B
1. YES
2. YES