NumPy Array Manipulation

Module 3: NumPy

NumPy: Array Manipulation

Numerical data does not always arrive in the shape a calculation needs. Imagine you play Fantasy Premier League (FPL) in a mini-league with your friends. These are the season statistics of the five midfielders in your squad, and you want them with one row per player:

Player	Goals	Assists	Clean Sheets	Defensive Contributions	Yellow Cards	Red Cards	Own Goals	Penalties Missed
Midfielder 1	15	10	12	2	3	0	0	1
Midfielder 2	6	7	10	14	6	0	0	0
Midfielder 3	2	3	11	22	9	1	1	0
Midfielder 4	9	11	8	1	2	0	0	0
Midfielder 5	4	12	9	5	4	0	0	1

However, the website you download them from exports the statistics stat by stat (all the goals, then all the assists…):

Statistic	Midfielder 1	Midfielder 2	Midfielder 3	Midfielder 4	Midfielder 5
Goals	15	6	2	9	4
Assists	10	7	3	11	12
Clean Sheets	12	10	11	8	9
Defensive Contributions	2	14	22	1	5
Yellow Cards	3	6	9	2	4
Red Cards	0	0	1	0	0
Own Goals	0	0	1	0	0
Penalties Missed	1	0	0	0	1

How do you turn the second table into the first? You need to reshape the array and transpose it.

And what happens if you sign a new player in the January transfer window? You need to add a new row to the array.

And if you want to compare your midfielders with those of a friend in your mini-league? You need to join two arrays.

The statistics can be grouped into attacking (goals and assists), defensive (clean sheets and defensive contributions), and disciplinary (yellow cards, red cards, own goals, and penalties missed). You may want to split the array into three arrays, one for each group.

In this lesson, we will learn how to manipulate arrays in NumPy, including reshaping, flattening, joining, and splitting arrays.

But first, as always, we need to import NumPy.

import numpy as np

Reshape Arrays

Any NumPy array has access to the reshape() method. This method creates a new NumPy array with the specified shape.

# We create an array of length 6
x = np.arange(0, 6, 1)
print(x)

[0 1 2 3 4 5]

# Reshape array to 2 x 3
reshaped = x.reshape(2, 3)
print(reshaped)
print(reshaped.shape)

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

The method reshape() creates a new array, it does not modify the original!

print(f"Original:\n{x}")
print(f"Reshaped:\n{reshaped}")

Original:
[0 1 2 3 4 5]
Reshaped:
[[0 1 2]
 [3 4 5]]

Exercise: Contact Sheet

A photographer prints a contact sheet of photos numbered 1 to 6, in 3 rows of 2. Create the following matrix of photo numbers without using the function np.array.

$\begin{pmatrix} 1 & 2\\ 3 & 4\\ 5 & 6 \end{pmatrix}$

Solution

arr = np.arange(1, 7).reshape(3, 2)
print(arr)

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

Exercise: Fill by Columns

The printer fills the sheet column by column instead. Create the following matrix without using the function np.array.

$\begin{pmatrix} 1 & 4\\ 2 & 5\\ 3 & 6 \end{pmatrix}$

Solution

arr = np.arange(1, 7).reshape(2, 3).T
print(arr)

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

Exercise: The Sudoku Box

A Sudoku grid has 81 cells. Number them from 0 to 80, row by row, in a $9 \times 9$ grid, without using the function np.array.

Then extract the central $3 \times 3$ box:

$\begin{pmatrix} 30 & 31 & 32\\ 39 & 40 & 41\\ 48 & 49 & 50\\ \end{pmatrix}$

Solution

grid = np.arange(81).reshape(9, 9)
print(grid)
print(grid[3:6, 3:6])

[[ 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 32 33 34 35]
 [36 37 38 39 40 41 42 43 44]
 [45 46 47 48 49 50 51 52 53]
 [54 55 56 57 58 59 60 61 62]
 [63 64 65 66 67 68 69 70 71]
 [72 73 74 75 76 77 78 79 80]]
[[30 31 32]
 [39 40 41]
 [48 49 50]]

When using reshape, you can add a single -1 to any dimension. When you do, NumPy decides the axis size to match the total number of elements in the array. For instance:

given an array with $N$ elements, calling reshape(a, -1) will return a shape (a, b) where b = N / a
given an array with $N$ elements, calling reshape(a, b, -1) will return a shape (a, b, c) where c = N / (a * b)

As the dimensions of reshape must always match exactly the number of elements in an array, using -1 will only work if the division produces an integer.

Reshape Is Not Transpose

The website exports the statistics stat by stat. Here are four statistics (goals, assists, clean sheets, and defensive contributions) of three players:

site_export = np.array([
    15, 6, 2,       # goals: player 1, player 2, player 3
    10, 7, 3,       # assists
    12, 10, 11,     # clean sheets
    2, 14, 22,      # defensive contributions
])

by_stat = site_export.reshape(4, 3)
print(by_stat)

[[15  6  2]
 [10  7  3]
 [12 10 11]
 [ 2 14 22]]

We want one row per player. A tempting mistake is to reshape to (3, 4):

print(site_export.reshape(3, 4))

[[15  6  2 10]
 [ 7  3 12 10]
 [11  2 14 22]]

The shape is right, but the values are wrong: reshape always fills row by row, so the first row mixes the goals and the assists. What we need is to swap rows and columns, the transpose .T:

by_player = by_stat.T
print(by_player)

[[15 10 12  2]
 [ 6  7 10 14]
 [ 2  3 11 22]]

Exercise: Spot the Difference

Without running any code, what is the first row of site_export.reshape(3, 4) and what is the first row of by_stat.T? Which one contains the statistics of player 1?

Solution

[15  6  2 10]
[15 10 12  2]

The first row of the reshape is the first four values of the export: the goals of the three players and the assists of player 1. The first row of the transpose is the first column of by_stat: the four statistics of player 1. Only the transpose contains the statistics of player 1.

Flatten Arrays

Another useful method is flatten(). It will create a one-dimensional array from a multi-dimensional array, effectively collapsing all the dimensions into a single dimension.

# We create a matrix 3x3
x = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
print(f"Original:\n{x}")

Original:
[[1 2 3]
 [4 5 6]
 [7 8 9]]

# Flatten array
flattened = x.flatten()
print(f"Flattened:\n{flattened}")

Flattened:
[1 2 3 4 5 6 7 8 9]

Joining Arrays

In data processing and analysis, you often need to combine arrays. This can be done in several ways depending on the structure and the desired outcome. Below, we’ll explore some common methods for joining arrays using NumPy.

Concatenation

Concatenation is the process of joining two or more arrays along an existing axis.

Concatenation on axis 0:

$\text{concat}_0(\mathbf{a}, \mathbf{b}) = [a_0, a_1, b_0, b_1]$

$_0(A, B) = \[\begin{pmatrix} a_{00} & a_{01} \\ a_{10} & a_{11} \\ b_{00} & b_{01} \\ b_{10} & b_{11} \\ \end{pmatrix}\]

# Create two 2-D arrays
array1 = np.array([[1, 2], [3, 4]])
array2 = np.array([[5, 6], [7, 8]])
print("Original shapes", array1.shape, "&", array2.shape)

# Concatenate along rows (axis=0)
concatenated_rows = np.concatenate((array1, array2), axis=0)
print("\nConcatenated along rows:\n", concatenated_rows)
print("\nConcatenated shape:", concatenated_rows.shape)

Original shapes (2, 2) & (2, 2)

Concatenated along rows:
 [[1 2]
 [3 4]
 [5 6]
 [7 8]]

Concatenated shape: (4, 2)

When you use np.concatenate(a, b, axis=0), it combines two arrays by stacking them vertically, one on top of the other. Imagine you have two bricks, and you place one brick directly on top of the other.

Now, let’s compare this to np.sum(a, axis=0). When you sum an array along axis=0, you are adding up the elements in each column. It’s like pressing down on the columns and “squeezing” the values together into a single row.

array1 = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
print("Original array:\n", array1)

result = np.sum(array1, axis=0)
print("We 'squeeze' the values down:\n  |  |  |\n  V  V  V\n", result)

Original array:
 [[1 2 3]
 [4 5 6]
 [7 8 9]]
We 'squeeze' the values down:
  |  |  |
  V  V  V
 [12 15 18]

Concatenation on axis 1:

$\text{concat}_1(\mathbf{a}, \mathbf{b}) = \text{ERROR}$

$_1(A, B) = \[\begin{pmatrix} a_{00} & a_{01} & b_{00} & b_{01} \\ a_{10} & a_{11} & b_{10} & b_{11} \\ \end{pmatrix}\]

# Create two 2-D arrays
array1 = np.array([[1, 2], [3, 4]])
array2 = np.array([[5, 6], [7, 8]])
print("Original shapes", array1.shape, "+", array2.shape)

# Concatenate along columns (axis=1)
concatenated_columns = np.concatenate((array1, array2), axis=1)
print("\nConcatenated along columns:\n", concatenated_columns)
print("\nConcatenated shape:", concatenated_columns.shape)

Original shapes (2, 2) + (2, 2)

Concatenated along columns:
 [[1 2 5 6]
 [3 4 7 8]]

Concatenated shape: (2, 4)

When you use np.concatenate(a, b, axis=1), it combines two arrays by stacking them horizontally, side to side. Imagine you have two bricks, and you place one brick at the left of the other.

Compare this to np.sum(a, axis=1). When you sum an array along axis=1, you are adding up the elements in each row. It’s like “squeezing” the values together into a single column.

array1 = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
print("Original array:\n", array1)

result = np.sum(array1, axis=1)
print("We 'squeeze' the values sideways:\n  ->\n  -> \n  ->\n", result)

Original array:
 [[1 2 3]
 [4 5 6]
 [7 8 9]]
We 'squeeze' the values sideways:
  ->
  -> 
  ->
 [ 6 15 24]

Pay attention to the following errors:

# You cannot concatenate on a non-existing axis
array1 = np.array([[1, 2], [3, 4]])
array2 = np.array([[5, 6], [7, 8]])
arr = np.concatenate((array1, array2), axis=2)

# The arrays must have the same number of dimensions
arr1 = np.array([[1, 2, 3], [4, 5, 6]])
arr2 = np.array([7, 8, 9])
arr = np.concatenate((arr1, arr2), axis=0)
print(arr)

Stacking

Stacking is another way to join arrays, but it can add an extra dimension to the result.

Horizontal stack: Stack arrays in sequence horizontally (column wise).

$\text{hstack}(\mathbf{a}, \mathbf{b}) = [a_0, a_1, b_0, b_1]$

$(A, B) = \[\begin{pmatrix} a_{00} & a_{01} & b_{00} & b_{01} \\ a_{10} & a_{11} & b_{10} & b_{11} \\ \end{pmatrix}\]

This is equivalent to concatenation along the second axis, except for 1-D arrays where it concatenates along the first axis.

# Create two 1-D arrays
array1 = np.array([1, 2, 3, 4])
array2 = np.array([5, 6, 7, 8])
print("Original shapes", array1.shape, "+", array2.shape)

# Horizontal stacking
hstacked = np.hstack((array1, array2))
print("\nHorizontally stacked:\n", hstacked)
print("\nHorizontally stacked shape:", hstacked.shape)

Original shapes (4,) + (4,)

Horizontally stacked:
 [1 2 3 4 5 6 7 8]

Horizontally stacked shape: (8,)

# Create two 2-D arrays
array1 = np.array([[1, 2], [3, 4]])
array2 = np.array([[5, 6], [7, 8]])
print("Original shapes", array1.shape, "+", array2.shape)

# Horizontal stacking
hstacked = np.hstack((array1, array2))
print("\nHorizontally stacked:\n", hstacked)
print("\nHorizontally stacked shape:", hstacked.shape)

Original shapes (2, 2) + (2, 2)

Horizontally stacked:
 [[1 2 5 6]
 [3 4 7 8]]

Horizontally stacked shape: (2, 4)

Vertical stack: Stack arrays in sequence vertically (row wise).

$\text{vstack}(\mathbf{a}, \mathbf{b}) = \begin{pmatrix} a_0 & a_1 \\ b_0 & b_1 \\ \end{pmatrix}$

$(A, B) = \[\begin{pmatrix} a_{00} & a_{01} \\ a_{10} & a_{11} \\ b_{00} & b_{01} \\ b_{10} & b_{11} \\ \end{pmatrix}\]

This is equivalent to concatenation along the first axis after 1-D arrays of shape (N,) have been reshaped to (1,N).

# Create two 1-D arrays
array1 = np.array([1, 2, 3, 4])
array2 = np.array([5, 6, 7, 8])
print("Original shapes", array1.shape, "+", array2.shape)

# Vertical stacking
vstacked = np.vstack((array1, array2))
print("\nVertically stacked:\n", vstacked)
print("\nVertically stacked shape:", vstacked.shape)

Original shapes (4,) + (4,)

Vertically stacked:
 [[1 2 3 4]
 [5 6 7 8]]

Vertically stacked shape: (2, 4)

# Create two 2-D arrays
array1 = np.array([[1, 2], [3, 4]])
array2 = np.array([[5, 6], [7, 8]])
print("Original shapes", array1.shape, "+", array2.shape)

# Horizontal stacking
vstacked = np.vstack((array1, array2))
print("\nVertically stacked:\n", vstacked)
print("\nVertically stacked shape:", vstacked.shape)

Original shapes (2, 2) + (2, 2)

Vertically stacked:
 [[1 2]
 [3 4]
 [5 6]
 [7 8]]

Vertically stacked shape: (4, 2)

# vstack is equilavent to ...

# Create two 1-D arrays
array1 = np.array([1, 2, 3, 4])
array2 = np.array([5, 6, 7, 8])
print("Original shapes", array1.shape, "+", array2.shape)

# Reshape the arrays to (1, 4)
array1 = array1.reshape(1, -1)
array2 = array2.reshape(1, -1)
print("Reshaped", array1.shape, "+", array2.shape)

# Concatenate along axis 1
concatenated_columns = np.concatenate((array1, array2), axis=0)
print("\nConcatenated along columns:\n", concatenated_columns)
print("\nConcatenated shape:", concatenated_columns.shape)

Original shapes (4,) + (4,)
Reshaped (1, 4) + (1, 4)

Concatenated along columns:
 [[1 2 3 4]
 [5 6 7 8]]

Concatenated shape: (2, 4)

Exercise: Rouché-Frobenius Theorem

A non-homogeneous system of linear equations $Ax=b$ with $n$ variables has a solution if and only if the rank of its coefficient matrix $A$ is equal to the rank of its augmented matrix $[A|b]$. If there are solutions, we can find that:

if $rank(A) = n$, the solution is unique,
if $rank(A) < n$, there are infinite solutions.

Write a function rouche_frobenius_theorem(A, b) that builds the augmented matrix $[A|b]$ with a stacking function and returns a message saying whether the system is inconsistent, has a unique solution, or has infinitely many solutions. Use np.linalg.matrix_rank() to compute ranks.

Test it with:

A = np.array([[2, 1, -1], [1, 3, 2], [1, -1, 2]])
b = np.array([8, 13, 3])

Solution

def rouche_frobenius_theorem(A, b):
    # Turn vector b with shape (N,) into a column with shape (N, 1).
    B = b.reshape(-1, 1)
    Ab = np.hstack((A, B))
    rank_A = np.linalg.matrix_rank(A)
    rank_Ab = np.linalg.matrix_rank(Ab)

    if rank_A == rank_Ab:
        if rank_A == A.shape[1]:
            return "The system is consistent and has a unique solution."
        return "The system is consistent and has infinitely many solutions."
    else:
        return "The system is inconsistent."

# Example usage
A = np.array([[2, 1, -1], [1, 3, 2], [1, -1, 2]])
b = np.array([8, 13, 3])

result = rouche_frobenius_theorem(A, b)
print(result)

The system is consistent and has a unique solution.

The vector b must be reshaped into a column before np.hstack() can place it next to A.

Stack

Column stacking is used to stack 1-D arrays as columns into a 2-D array.

$\text{column\_stack}(\mathbf{a}, \mathbf{b}) = \begin{pmatrix} a_0 & b_0 \\ a_1 & b_1 \\ \end{pmatrix}$

Take a sequence of 1-D arrays and stack them as columns to make a single 2-D array. 2-D arrays are stacked as-is, just like with hstack. 1-D arrays are turned into 2-D columns first.

# Create two 1-D arrays
array1 = np.array([1, 2, 3])
array2 = np.array([4, 5, 6])

# Column stack
column_stacked = np.column_stack((array1, array2))
print("Column stacked:\n", column_stacked)

Column stacked:
 [[1 4]
 [2 5]
 [3 6]]

Row stacking is used to stack 1-D arrays as rows into a 2-D array.

$\text{row\_stack}(\mathbf{a}, \mathbf{b}) = \begin{pmatrix} a_0 & a_1 \\ b_0 & b_1 \\ \end{pmatrix}$

This is equivalent to concatenation along the first axis after 1-D arrays of shape (N,) have been reshaped to (1,N). In older code you will find it as np.row_stack(), an alias of np.vstack() that is deprecated since NumPy 2.0. Use np.vstack() instead.

# Create two 1-D arrays
array1 = np.array([1, 2, 3])
array2 = np.array([4, 5, 6])

# Row stack
row_stacked = np.vstack((array1, array2))
print("Row stacked:\n", row_stacked)

Row stacked:
 [[1 2 3]
 [4 5 6]]

Back to the league. This is the first table of the introduction as a NumPy array: one row per player, and one column for each of the eight statistics.

# Season statistics of our midfielders: one row per player, one column per statistic
# goals, assists, clean sheets, defensive contributions,
# yellow cards, red cards, own goals, penalties missed
stats = np.array([
    [15, 10, 12,  2, 3, 0, 0, 1],
    [ 6,  7, 10, 14, 6, 0, 0, 0],
    [ 2,  3, 11, 22, 9, 1, 1, 0],
    [ 9, 11,  8,  1, 2, 0, 0, 0],
    [ 4, 12,  9,  5, 4, 0, 0, 1],
])

Exercise: The January Transfer

You sign a new midfielder in the January transfer window. Add his statistics to stats as a new row, and print the shape of the result.

new_player = np.array([5, 4, 6, 3, 2, 0, 0, 0])

Solution

new_player = np.array([5, 4, 6, 3, 2, 0, 0, 0])

with_new = np.vstack((stats, new_player))
print(with_new.shape)

(6, 8)

new_player has shape (8,), so np.vstack() treats it as a row of shape (1, 8) and places it below the 5 rows of stats: the result has shape (6, 8).

Splitting Arrays

Splitting arrays is the process of dividing an array into multiple sub-arrays. This can be useful for various data manipulation tasks. Below, we’ll explore some common methods for splitting arrays using NumPy.

Split

The np.split() function splits an array into multiple sub-arrays along a specified axis.

# Create an array
arr = np.array([1, 2, 3, 4, 5, 6])

# Split the array into 3 equal parts
split_array = np.split(arr, 3)

for idx in range(len(split_array)):
    print(f"Split {idx + 1}:\n{split_array[idx]}")

Split 1:
[1 2]
Split 2:
[3 4]
Split 3:
[5 6]

Array Split

The np.array_split() function allows you to split an array into unequal parts if needed.

For an array of length l that should be split into n sections, it returns l % n sub-arrays of size l//n + 1 and the rest of size l//n.

# Split the array into 4 parts (unequal)
split_array_unequal = np.array_split(arr, 4)

for idx in range(len(split_array_unequal)):
    print(f"Split {idx + 1}:\n{split_array_unequal[idx]}")

Split 1:
[1 2]
Split 2:
[3 4]
Split 3:
[5]
Split 4:
[6]

Horizontal Split

The np.hsplit() function splits an array horizontally (column-wise).

# Create a 2D array
array_2d = np.array([[1, 2, 3], [4, 5, 6]])

# Split the 2D array into 3 columns
hsplit_array = np.hsplit(array_2d, 3)

for idx in range(len(hsplit_array)):
    print(f"Split {idx + 1}:\n{hsplit_array[idx]}")

Split 1:
[[1]
 [4]]
Split 2:
[[2]
 [5]]
Split 3:
[[3]
 [6]]

Vertical Split

The np.vsplit() function splits an array vertically (row-wise).

# Split the 2D array into 2 rows
vsplit_array = np.vsplit(array_2d, 2)

for idx in range(len(vsplit_array)):
    print(f"Split {idx + 1}:\n{vsplit_array[idx]}")

Split 1:
[[1 2 3]]
Split 2:
[[4 5 6]]

Splitting at Given Positions

Instead of a number of equal parts, all the split functions also accept a list of positions where to cut. The three groups of statistics have different sizes: attack (goals and assists), defence (clean sheets and defensive contributions), and discipline (yellow cards, red cards, own goals, and penalties missed). We cut the columns before index 2 and before index 4:

attack, defence, discipline = np.hsplit(stats, [2, 4])

print(attack.shape, defence.shape, discipline.shape)

(5, 2) (5, 2) (5, 4)

In FPL, a midfielder gets 5 points for each goal, 3 for each assist, 1 for each clean sheet, and 2 for each match with defensive contributions. He loses 1 point for each yellow card, 3 for each red card, 2 for each own goal, and 2 for each penalty missed. We store the points in a one-dimensional array, in the same order as the columns of stats:

# FPL points of a midfielder for each statistic
points = np.array([5, 3, 1, 2, -1, -3, -2, -2])

Exercise: Split the Points

Split points into three arrays, one for each group of statistics: attack (the first two values), defence (the next two), and discipline (the last four). Then check your split by computing the sum of each group: one goal and one assist are worth 8 points, one clean sheet and one defensive contribution are worth 3 points, and one of each disciplinary event costs 8 points.

Solution

attack_p, defence_p, discipline_p = np.split(points, [2, 4])

print(np.sum(attack_p), np.sum(defence_p), np.sum(discipline_p))

8 3 -8

points is one-dimensional, so we use np.split() (or np.hsplit(), which behaves the same for 1D arrays).

Challenge: One Array for the Whole League

A friend in your mini-league sent you the statistics of the 5 midfielders in their squad stat by stat, as one long array: the 5 goals, then the 5 assists, and so on. They follow the same order of statistics as our team.

rival_export = np.array([
    11, 3, 7, 1, 5,       # goals
    6, 5, 9, 2, 8,        # assists
    10, 12, 7, 13, 9,     # clean sheets
    3, 16, 2, 25, 6,      # defensive contributions
    2, 5, 1, 8, 3,        # yellow cards
    0, 0, 0, 1, 0,        # red cards
    0, 1, 0, 0, 0,        # own goals
    1, 0, 0, 0, 1,        # penalties missed
])

Without writing any loop:

Turn rival_export into an array with one row per player and one column per statistic.
Join it with our stats into a single array, league, with our players first. Print its shape.
Split league into the three groups of statistics, and compute the fantasy points of each player in each group (use points).
Compute the total fantasy points of every player.
Build a report with two columns: the player number (1 to 10) and the total points.

Possible solution

# 1. One row per statistic first, then transpose
rival = rival_export.reshape(8, 5).T
print(rival)

# 2. Join along the rows
league = np.concatenate((stats, rival), axis=0)  # or np.vstack()
print(league.shape)

# 3. Split the columns at the start of each group
attack, defence, discipline = np.hsplit(league, [2, 4])
attack_p, defence_p, discipline_p = np.split(points, [2, 4])

print(np.sum(attack * attack_p, axis=1))
print(np.sum(defence * defence_p, axis=1))
print(np.sum(discipline * discipline_p, axis=1))

# 4. Total fantasy points
total = league @ points
print(total)

# 5. Report
report = np.column_stack((np.arange(1, 11), total))
print(report)

[[11  6 10  3  2  0  0  1]
 [ 3  5 12 16  5  0  1  0]
 [ 7  9  7  2  1  0  0  0]
 [ 1  2 13 25  8  1  0  0]
 [ 5  8  9  6  3  0  0  1]]
(10, 8)
[105  51  19  78  56  73  30  62  11  49]
[16 38 55 10 19 16 44 11 63 21]
[ -5  -6 -14  -2  -6  -4  -7  -1 -11  -5]
[116  83  60  86  69  85  67  72  63  65]
[[  1 116]
 [  2  83]
 [  3  60]
 [  4  86]
 [  5  69]
 [  6  85]
 [  7  67]
 [  8  72]
 [  9  63]
 [ 10  65]]

Note that rival_export.reshape(5, 8) would also give the right shape, but with the statistics mixed up.

Exam-Style Practice

The midterm and the final exam have the same kinds of exercises you have been solving: a debugging exercise and an open problem. Try these two without AI, as you would in the exam, and open the solution only when you are done. Both come from the course of Manuel Morante.

Debugging Exercise: Padding a Matrix

The following function should add a border of a given number around a matrix, like the function np.pad(). It adds one row above, one row below, one column on the left, and one column on the right. This can be useful when every element needs neighbors in all directions.

However, it raises this error:

ValueError: all the input arrays must have same number of dimensions, but the array at index 0 has 1 dimension(s) and the array at index 1 has 2 dimension(s)

Explain the error and fix the function. Does your fix work for a matrix that is not square?

def manual_pad(A, number):
    # Create rows and columns to add around the current array
    rows = np.ones(A.shape[0]) * number
    columns = np.ones(A.shape[1]) * number

    # Stack vertically and horizontally the new rows and columns
    C = np.vstack((rows, A, rows))
    D = np.hstack((columns, C, columns))

    return D

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

print(manual_pad(A, 1))

Solution

There are three problems, and each one is about shapes:

np.vstack() accepts the 1D rows because it treats a vector as a row, so the vertical stack works.
np.hstack() needs arrays with the same number of dimensions. columns is 1D, and C is 2D. This is the error. A column must have the shape (n, 1), not (n,).
The new columns must be as long as the array after adding the rows, which is C.shape[0], not A.shape[1]. Also, a new row must have A.shape[1] values, not A.shape[0]. The code only works by chance when the matrix is square, and a rectangular matrix would fail with a different error.

def manual_pad(A, number):
    # One row with as many values as columns, shape (1, columns)
    rows = np.full((1, A.shape[1]), number)
    C = np.vstack((rows, A, rows))

    # One column with as many values as the new number of rows, shape (rows, 1)
    columns = np.full((C.shape[0], 1), number)
    return np.hstack((columns, C, columns))


A = np.array([[1, 0, -1],
              [-2, 1, 2],
              [0, -3, 1]])
R = np.array([[1, 2, 3, 4],
              [5, 6, 7, 8]])   # not square

print(manual_pad(A, 1))
print(manual_pad(R, 0))

[[ 1  1  1  1  1]
 [ 1  1  0 -1  1]
 [ 1 -2  1  2  1]
 [ 1  0 -3  1  1]
 [ 1  1  1  1  1]]
[[0 0 0 0 0 0]
 [0 1 2 3 4 0]
 [0 5 6 7 8 0]
 [0 0 0 0 0 0]]

NumPy already has this function. Check that your version agrees with it:

print(np.array_equal(manual_pad(R, 7), np.pad(R, 1, constant_values=7)))

True

Knowing that np.pad() exists is what matters in the exam, because you will find it in the cheat sheet. Understanding the shapes is what lets you use it, or write it yourself when it is not there.

Open Problem: Vectorized Evaluation of Mathematical Functions

Your mathematics teacher asks you to apply several functions to a range of numbers. The challenge is to use a single array for the whole range of input values. Do not create a separate input array for each function.

Create a 1D array with 21 equally spaced values from 0.5 to 4, including both endpoints. Its length is divisible by 3.
Split it into three equal sub-arrays. Evaluate the following functions on the corresponding sub-array, with vectorized operations and no loops:
- $f_1(x) = x^3 - 3x + 2$ on the first sub-array.
- $f_2(x) = \ln(x) + \sqrt{x}$ on the second sub-array.
- $f_3(x) = e^{x/2}$ on the third sub-array.
For each function, find the positions where its result is greater than 2. A position is the index inside its own sub-array, from 0 to 6.
Build a 2D array with three columns, one for each function, in this order. Each column has the positions of its function. The columns may have different lengths, so fill the shorter ones with np.nan.

Solution

x = np.linspace(0.5, 4, 21)
x1, x2, x3 = np.split(x, 3)

f1 = x1**3 - 3*x1 + 2
f2 = np.log(x2) + np.sqrt(x2)
f3 = np.exp(x3 / 2)

# One array of positions for each function
positions = [np.where(f > 2)[0] for f in (f1, f2, f3)]
print(positions)

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

The columns have different lengths, so we first create the final array full of np.nan and then fill each column with slicing. The number of rows is the length of the longest column:

n_rows = max(len(p) for p in positions)

result = np.full((n_rows, 3), np.nan)
for column, p in enumerate(positions):
    result[:len(p), column] = p

print(result)

[[nan  1.  0.]
 [nan  2.  1.]
 [nan  3.  2.]
 [nan  4.  3.]
 [nan  5.  4.]
 [nan  6.  5.]
 [nan nan  6.]]

On this range, $f_1$ never exceeds 2 (its largest value is about 1.07), so np.where returns an empty array and the first column is all np.nan. A correct solution must not fail in that case. The loop above only assembles the result: the functions themselves are applied without loops.

Exam-Style Practice

The midterm and the final exam have the same kinds of exercises you have been solving: a debugging exercise and an open problem. Try each one on your own, without AI, as you would in the exam, before opening its collapsed solution.

Debugging: Padding a Matrix

Adapted from the exam-style exercises of Manuel Morante.

The following function should add a border of a given number around a matrix, like the function np.pad(). It adds one row above, one row below, one column on the left, and one column on the right. This is useful when every element needs neighbors in all directions.

However, it raises this error:

ValueError: all the input arrays must have same number of dimensions, but the array at index 0 has 1 dimension(s) and the array at index 1 has 2 dimension(s)

Explain the error and fix the function. Does your fix work for a matrix that is not square?

def manual_pad(A, number):
    # Create rows and columns to add around the current array
    rows = np.ones(A.shape[0]) * number
    columns = np.ones(A.shape[1]) * number

    # Stack vertically and horizontally the new rows and columns
    C = np.vstack((rows, A, rows))
    D = np.hstack((columns, C, columns))

    return D

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

print(manual_pad(A, 1))

Solution

There are three problems, and each one is about shapes:

np.vstack() accepts the 1D rows because it treats a vector as a row, so the vertical stack works.
np.hstack() needs arrays with the same number of dimensions. columns is 1D, and C is 2D. This is the error. A column must have the shape (n, 1), not (n,).
The sizes are wrong as well. A new row must have A.shape[1] values (one per column), not A.shape[0], and a new column must be as long as the array after adding the rows, which is C.shape[0]. The original code only has matching sizes for a square matrix, so a rectangular one would still fail after fixing the dimensions.

def manual_pad(A, number):
    # One row with as many values as columns, shape (1, columns)
    rows = np.full((1, A.shape[1]), number)
    C = np.vstack((rows, A, rows))

    # One column with as many values as the new number of rows, shape (rows, 1)
    columns = np.full((C.shape[0], 1), number)
    return np.hstack((columns, C, columns))


A = np.array([[1, 0, -1],
              [-2, 1, 2],
              [0, -3, 1]])
R = np.array([[1, 2, 3, 4],
              [5, 6, 7, 8]])   # not square

print(manual_pad(A, 1))
print(manual_pad(R, 0))

[[ 1  1  1  1  1]
 [ 1  1  0 -1  1]
 [ 1 -2  1  2  1]
 [ 1  0 -3  1  1]
 [ 1  1  1  1  1]]
[[0 0 0 0 0 0]
 [0 1 2 3 4 0]
 [0 5 6 7 8 0]
 [0 0 0 0 0 0]]

NumPy already has this function. Check that your version agrees with it:

print(np.array_equal(manual_pad(R, 7), np.pad(R, 1, constant_values=7)))

True

Knowing that np.pad() exists is what matters in the exam, because you will find it in the cheat sheet. Understanding the shapes is what lets you use it, or write it yourself when it is not there.

Open Problem: Vectorized Evaluation of Mathematical Functions

Adapted from the exam-style exercises of Manuel Morante.

Create a 1D array with 21 equally spaced values from 0.5 to 4, including both endpoints. Its length is divisible by 3.
Split it into three equal sub-arrays. Evaluate the following functions on the corresponding sub-array, with vectorized operations and no loops to apply them:
- $f_1(x) = x^3 - 3x + 2$ on the first sub-array.
- $f_2(x) = \ln(x) + \sqrt{x}$ on the second sub-array.
- $f_3(x) = e^{x/2}$ on the third sub-array.
For each function, find the positions where its result is greater than 2. A position is the index inside its own sub-array, from 0 to 6.
Build a 2D array with three columns, one for each function, in this order. Each column has the positions of its function. The columns may have different lengths, so fill the shorter ones with np.nan.

Expected result:

[[nan  1.  0.]
 [nan  2.  1.]
 [nan  3.  2.]
 [nan  4.  3.]
 [nan  5.  4.]
 [nan  6.  5.]
 [nan nan  6.]]

Solution

x = np.linspace(0.5, 4, 21)
x1, x2, x3 = np.split(x, 3)

f1 = x1**3 - 3*x1 + 2
f2 = np.log(x2) + np.sqrt(x2)
f3 = np.exp(x3 / 2)

# One array of positions for each function
positions = [np.where(f > 2)[0] for f in (f1, f2, f3)]
print(positions)

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

The columns have different lengths, so we first create the final array full of np.nan and then fill each column with slicing. The number of rows is the length of the longest column:

n_rows = max(len(p) for p in positions)

result = np.full((n_rows, 3), np.nan)
for column, p in enumerate(positions):
    result[:len(p), column] = p

print(result)

[[nan  1.  0.]
 [nan  2.  1.]
 [nan  3.  2.]
 [nan  4.  3.]
 [nan  5.  4.]
 [nan  6.  5.]
 [nan nan  6.]]

Summary

You do not need to memorize every joining and splitting function. Focus on the shape of each input, the shape you want as output, and the axis along which the operation acts.

reshape() changes the arrangement of values but not their count.
flatten() returns a one-dimensional copy.
concatenate() joins arrays along an existing axis.
stack(), vstack(), and hstack() provide common arrangements.
split() requires equal pieces; array_split() also permits unequal pieces.

Before running any joining operation, predict the output shape. Then explain why all non-joining dimensions must match.

Homework and Review

Work through the homework exercises available here.

Homework

--- title: "NumPy Array Manipulation" subtitle: "Module 3: NumPy" format: html --- # NumPy: Array Manipulation Numerical data does not always arrive in the shape a calculation needs. Imagine you play **Fantasy Premier League** (FPL) in a mini-league with your friends. These are the season statistics of the five midfielders in your squad, and you want them with one row per player: | Player | Goals | Assists | Clean Sheets | Defensive Contributions | Yellow Cards | Red Cards | Own Goals | Penalties Missed | |--------|-------|---------|--------------|-------------------------|--------------|-----------|-----------|------------------| | Midfielder 1 | 15 | 10 | 12 | 2 | 3 | 0 | 0 | 1 | | Midfielder 2 | 6 | 7 | 10 | 14 | 6 | 0 | 0 | 0 | | Midfielder 3 | 2 | 3 | 11 | 22 | 9 | 1 | 1 | 0 | | Midfielder 4 | 9 | 11 | 8 | 1 | 2 | 0 | 0 | 0 | | Midfielder 5 | 4 | 12 | 9 | 5 | 4 | 0 | 0 | 1 | However, the website you download them from exports the statistics **stat by stat** (all the goals, then all the assists...): | Statistic | Midfielder 1 | Midfielder 2 | Midfielder 3 | Midfielder 4 | Midfielder 5 | |-----------|--------------|--------------|--------------|--------------|--------------| | Goals | 15 | 6 | 2 | 9 | 4 | | Assists | 10 | 7 | 3 | 11 | 12 | | Clean Sheets | 12 | 10 | 11 | 8 | 9 | | Defensive Contributions | 2 | 14 | 22 | 1 | 5 | | Yellow Cards | 3 | 6 | 9 | 2 | 4 | | Red Cards | 0 | 0 | 1 | 0 | 0 | | Own Goals | 0 | 0 | 1 | 0 | 0 | | Penalties Missed | 1 | 0 | 0 | 0 | 1 | How do you turn the second table into the first? You need to **reshape** the array and **transpose** it. And what happens if you sign a new player in the January transfer window? You need to **add a new row** to the array. And if you want to compare your midfielders with those of a friend in your mini-league? You need to **join two arrays**. The statistics can be grouped into attacking (goals and assists), defensive (clean sheets and defensive contributions), and disciplinary (yellow cards, red cards, own goals, and penalties missed). You may want to **split** the array into three arrays, one for each group. In this lesson, we will learn how to manipulate arrays in NumPy, including reshaping, flattening, joining, and splitting arrays. But first, as always, we need to import NumPy. ```{python} import numpy as np ``` ## Reshape Arrays Any NumPy array has access to the `reshape()` method. This method creates a new NumPy array with the specified shape. ```{python} # We create an array of length 6 x = np.arange(0, 6, 1) print(x) ``` ```{python} # Reshape array to 2 x 3 reshaped = x.reshape(2, 3) print(reshaped) print(reshaped.shape) ``` The method `reshape()` creates a new array, it does not modify the original! ```{python} print(f"Original:\n{x}") print(f"Reshaped:\n{reshaped}") ``` ::: {.callout-note} ## Exercise: Contact Sheet A photographer prints a contact sheet of photos numbered 1 to 6, in 3 rows of 2. Create the following matrix of photo numbers without using the function `np.array`. $\begin{pmatrix} 1 & 2\\ 3 & 4\\ 5 & 6 \end{pmatrix}$ ::: ::: {.callout-tip collapse="true"} ## Solution ```{python} arr = np.arange(1, 7).reshape(3, 2) print(arr) ``` ::: ::: {.callout-note} ## Exercise: Fill by Columns The printer fills the sheet **column by column** instead. Create the following matrix without using the function `np.array`. $\begin{pmatrix} 1 & 4\\ 2 & 5\\ 3 & 6 \end{pmatrix}$ ::: ::: {.callout-tip collapse="true"} ## Solution ```{python} arr = np.arange(1, 7).reshape(2, 3).T print(arr) ``` ::: ::: {.callout-note} ## Exercise: The Sudoku Box A Sudoku grid has 81 cells. Number them from 0 to 80, row by row, in a $9 \times 9$ grid, without using the function `np.array`. Then extract the **central** $3 \times 3$ box: $\begin{pmatrix} 30 & 31 & 32\\ 39 & 40 & 41\\ 48 & 49 & 50\\ \end{pmatrix}$ ::: ::: {.callout-tip collapse="true"} ## Solution ```{python} grid = np.arange(81).reshape(9, 9) print(grid) print(grid[3:6, 3:6]) ``` ::: When using reshape, you can add a single -1 to any dimension. When you do, NumPy decides the axis size to match the total number of elements in the array. For instance: - given an array with $N$ elements, calling `reshape(a, -1)` will return a shape `(a, b)` where `b = N / a` - given an array with $N$ elements, calling `reshape(a, b, -1)` will return a shape `(a, b, c)` where `c = N / (a * b)` As the dimensions of `reshape` must always match exactly the number of elements in an array, using -1 will only work if the division produces an integer. ### Reshape Is Not Transpose The website exports the statistics **stat by stat**. Here are four statistics (goals, assists, clean sheets, and defensive contributions) of three players: ```{python} site_export = np.array([ 15, 6, 2, # goals: player 1, player 2, player 3 10, 7, 3, # assists 12, 10, 11, # clean sheets 2, 14, 22, # defensive contributions ]) by_stat = site_export.reshape(4, 3) print(by_stat) ``` We want one row per player. A tempting mistake is to reshape to `(3, 4)`: ```{python} print(site_export.reshape(3, 4)) ``` The shape is right, but the values are wrong: reshape always fills row by row, so the first row mixes the goals and the assists. What we need is to swap rows and columns, the **transpose** `.T`: ```{python} by_player = by_stat.T print(by_player) ``` ::: {.callout-note} ## Exercise: Spot the Difference Without running any code, what is the first row of `site_export.reshape(3, 4)` and what is the first row of `by_stat.T`? Which one contains the statistics of player 1? ::: ::: {.callout-tip collapse="true"} ## Solution ```text [15 6 2 10] [15 10 12 2] ``` The first row of the reshape is the first four values of the export: the goals of the three players and the assists of player 1. The first row of the transpose is the first **column** of `by_stat`: the four statistics of player 1. Only the transpose contains the statistics of player 1. ::: --- ## Flatten Arrays Another useful method is `flatten()`. It will create a one-dimensional array from a multi-dimensional array, effectively collapsing all the dimensions into a single dimension. ```{python} # We create a matrix 3x3 x = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]]) print(f"Original:\n{x}") ``` ```{python} # Flatten array flattened = x.flatten() print(f"Flattened:\n{flattened}") ``` --- ## Joining Arrays In data processing and analysis, you often need to combine arrays. This can be done in several ways depending on the structure and the desired outcome. Below, we'll explore some common methods for joining arrays using NumPy. ### Concatenation Concatenation is the process of joining two or more arrays along an existing axis. **Concatenation on axis 0**: $\text{concat}_0(\mathbf{a}, \mathbf{b}) = [a_0, a_1, b_0, b_1]$ $\text{concat}_0(A, B) = \begin{pmatrix} a_{00} & a_{01} \\ a_{10} & a_{11} \\ b_{00} & b_{01} \\ b_{10} & b_{11} \\ \end{pmatrix} $ ```{python} # Create two 2-D arrays array1 = np.array([[1, 2], [3, 4]]) array2 = np.array([[5, 6], [7, 8]]) print("Original shapes", array1.shape, "&", array2.shape) # Concatenate along rows (axis=0) concatenated_rows = np.concatenate((array1, array2), axis=0) print("\nConcatenated along rows:\n", concatenated_rows) print("\nConcatenated shape:", concatenated_rows.shape) ``` When you use `np.concatenate(a, b, axis=0)`, it combines two arrays by stacking them vertically, one on top of the other. Imagine you have two bricks, and you place one brick directly on top of the other. Now, let's compare this to `np.sum(a, axis=0)`. When you sum an array along `axis=0`, you are adding up the elements in each column. It's like pressing down on the columns and "squeezing" the values together into a single row. ```{python} array1 = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]]) print("Original array:\n", array1) result = np.sum(array1, axis=0) print("We 'squeeze' the values down:\n | | |\n V V V\n", result) ``` **Concatenation on axis 1**: $\text{concat}_1(\mathbf{a}, \mathbf{b}) = \text{ERROR}$ $\text{concat}_1(A, B) = \begin{pmatrix} a_{00} & a_{01} & b_{00} & b_{01} \\ a_{10} & a_{11} & b_{10} & b_{11} \\ \end{pmatrix} $ ```{python} # Create two 2-D arrays array1 = np.array([[1, 2], [3, 4]]) array2 = np.array([[5, 6], [7, 8]]) print("Original shapes", array1.shape, "+", array2.shape) # Concatenate along columns (axis=1) concatenated_columns = np.concatenate((array1, array2), axis=1) print("\nConcatenated along columns:\n", concatenated_columns) print("\nConcatenated shape:", concatenated_columns.shape) ``` When you use `np.concatenate(a, b, axis=1)`, it combines two arrays by stacking them horizontally, side to side. Imagine you have two bricks, and you place one brick at the left of the other. Compare this to `np.sum(a, axis=1)`. When you sum an array along `axis=1`, you are adding up the elements in each row. It's like "squeezing" the values together into a single column. ```{python} array1 = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]]) print("Original array:\n", array1) result = np.sum(array1, axis=1) print("We 'squeeze' the values sideways:\n ->\n -> \n ->\n", result) ``` Pay attention to the following errors: ```{python} #| eval: false # You cannot concatenate on a non-existing axis array1 = np.array([[1, 2], [3, 4]]) array2 = np.array([[5, 6], [7, 8]]) arr = np.concatenate((array1, array2), axis=2) ``` ```{python} #| eval: false # The arrays must have the same number of dimensions arr1 = np.array([[1, 2, 3], [4, 5, 6]]) arr2 = np.array([7, 8, 9]) arr = np.concatenate((arr1, arr2), axis=0) print(arr) ``` ### Stacking Stacking is another way to join arrays, but it can add an extra dimension to the result. **Horizontal stack:** Stack arrays in sequence horizontally (column wise). $\text{hstack}(\mathbf{a}, \mathbf{b}) = [a_0, a_1, b_0, b_1]$ $\text{hstack}(A, B) = \begin{pmatrix} a_{00} & a_{01} & b_{00} & b_{01} \\ a_{10} & a_{11} & b_{10} & b_{11} \\ \end{pmatrix} $ This is equivalent to concatenation along the second axis, except for 1-D arrays where it concatenates along the first axis. ```{python} # Create two 1-D arrays array1 = np.array([1, 2, 3, 4]) array2 = np.array([5, 6, 7, 8]) print("Original shapes", array1.shape, "+", array2.shape) # Horizontal stacking hstacked = np.hstack((array1, array2)) print("\nHorizontally stacked:\n", hstacked) print("\nHorizontally stacked shape:", hstacked.shape) ``` ```{python} # Create two 2-D arrays array1 = np.array([[1, 2], [3, 4]]) array2 = np.array([[5, 6], [7, 8]]) print("Original shapes", array1.shape, "+", array2.shape) # Horizontal stacking hstacked = np.hstack((array1, array2)) print("\nHorizontally stacked:\n", hstacked) print("\nHorizontally stacked shape:", hstacked.shape) ``` **Vertical stack:** Stack arrays in sequence vertically (row wise). $\text{vstack}(\mathbf{a}, \mathbf{b}) = \begin{pmatrix} a_0 & a_1 \\ b_0 & b_1 \\ \end{pmatrix}$ $\text{vstack}(A, B) = \begin{pmatrix} a_{00} & a_{01} \\ a_{10} & a_{11} \\ b_{00} & b_{01} \\ b_{10} & b_{11} \\ \end{pmatrix} $ This is equivalent to concatenation along the first axis after 1-D arrays of shape _(N,)_ have been reshaped to _(1,N)_. ```{python} # Create two 1-D arrays array1 = np.array([1, 2, 3, 4]) array2 = np.array([5, 6, 7, 8]) print("Original shapes", array1.shape, "+", array2.shape) # Vertical stacking vstacked = np.vstack((array1, array2)) print("\nVertically stacked:\n", vstacked) print("\nVertically stacked shape:", vstacked.shape) ``` ```{python} # Create two 2-D arrays array1 = np.array([[1, 2], [3, 4]]) array2 = np.array([[5, 6], [7, 8]]) print("Original shapes", array1.shape, "+", array2.shape) # Horizontal stacking vstacked = np.vstack((array1, array2)) print("\nVertically stacked:\n", vstacked) print("\nVertically stacked shape:", vstacked.shape) ``` ```{python} # vstack is equilavent to ... # Create two 1-D arrays array1 = np.array([1, 2, 3, 4]) array2 = np.array([5, 6, 7, 8]) print("Original shapes", array1.shape, "+", array2.shape) # Reshape the arrays to (1, 4) array1 = array1.reshape(1, -1) array2 = array2.reshape(1, -1) print("Reshaped", array1.shape, "+", array2.shape) # Concatenate along axis 1 concatenated_columns = np.concatenate((array1, array2), axis=0) print("\nConcatenated along columns:\n", concatenated_columns) print("\nConcatenated shape:", concatenated_columns.shape) ``` ::: {.callout-note} ## Exercise: Rouché-Frobenius Theorem A non-homogeneous system of linear equations $Ax=b$ with $n$ variables has a solution if and only if the rank of its coefficient matrix $A$ is equal to the rank of its augmented matrix $[A|b]$. If there are solutions, we can find that: - if $rank(A) = n$, the solution is unique, - if $rank(A) < n$, there are infinite solutions. Write a function `rouche_frobenius_theorem(A, b)` that builds the augmented matrix $[A|b]$ with a stacking function and returns a message saying whether the system is inconsistent, has a unique solution, or has infinitely many solutions. Use `np.linalg.matrix_rank()` to compute ranks. Test it with: ```python A = np.array([[2, 1, -1], [1, 3, 2], [1, -1, 2]]) b = np.array([8, 13, 3]) ``` ::: ::: {.callout-tip collapse="true"} ## Solution ```{python} def rouche_frobenius_theorem(A, b): # Turn vector b with shape (N,) into a column with shape (N, 1). B = b.reshape(-1, 1) Ab = np.hstack((A, B)) rank_A = np.linalg.matrix_rank(A) rank_Ab = np.linalg.matrix_rank(Ab) if rank_A == rank_Ab: if rank_A == A.shape[1]: return "The system is consistent and has a unique solution." return "The system is consistent and has infinitely many solutions." else: return "The system is inconsistent." # Example usage A = np.array([[2, 1, -1], [1, 3, 2], [1, -1, 2]]) b = np.array([8, 13, 3]) result = rouche_frobenius_theorem(A, b) print(result) ``` The vector `b` must be reshaped into a column before `np.hstack()` can place it next to `A`. ::: ### Stack **Column stacking** is used to stack 1-D arrays as columns into a 2-D array. $\text{column\_stack}(\mathbf{a}, \mathbf{b}) = \begin{pmatrix} a_0 & b_0 \\ a_1 & b_1 \\ \end{pmatrix}$ Take a sequence of 1-D arrays and stack them as columns to make a single 2-D array. 2-D arrays are stacked as-is, just like with hstack. 1-D arrays are turned into 2-D columns first. ```{python} # Create two 1-D arrays array1 = np.array([1, 2, 3]) array2 = np.array([4, 5, 6]) # Column stack column_stacked = np.column_stack((array1, array2)) print("Column stacked:\n", column_stacked) ``` **Row stacking** is used to stack 1-D arrays as rows into a 2-D array. $\text{row\_stack}(\mathbf{a}, \mathbf{b}) = \begin{pmatrix} a_0 & a_1 \\ b_0 & b_1 \\ \end{pmatrix}$ This is equivalent to concatenation along the first axis after 1-D arrays of shape _(N,)_ have been reshaped to _(1,N)_. In older code you will find it as `np.row_stack()`, an alias of `np.vstack()` that is deprecated since NumPy 2.0. Use `np.vstack()` instead. ```{python} # Create two 1-D arrays array1 = np.array([1, 2, 3]) array2 = np.array([4, 5, 6]) # Row stack row_stacked = np.vstack((array1, array2)) print("Row stacked:\n", row_stacked) ``` Back to the league. This is the first table of the introduction as a NumPy array: one row per player, and one column for each of the eight statistics. ```{python} # Season statistics of our midfielders: one row per player, one column per statistic # goals, assists, clean sheets, defensive contributions, # yellow cards, red cards, own goals, penalties missed stats = np.array([ [15, 10, 12, 2, 3, 0, 0, 1], [ 6, 7, 10, 14, 6, 0, 0, 0], [ 2, 3, 11, 22, 9, 1, 1, 0], [ 9, 11, 8, 1, 2, 0, 0, 0], [ 4, 12, 9, 5, 4, 0, 0, 1], ]) ``` ::: {.callout-note} ## Exercise: The January Transfer You sign a new midfielder in the January transfer window. Add his statistics to `stats` as a new row, and print the shape of the result. ```python new_player = np.array([5, 4, 6, 3, 2, 0, 0, 0]) ``` ::: ::: {.callout-tip collapse="true"} ## Solution ```{python} new_player = np.array([5, 4, 6, 3, 2, 0, 0, 0]) with_new = np.vstack((stats, new_player)) print(with_new.shape) ``` `new_player` has shape `(8,)`, so `np.vstack()` treats it as a row of shape `(1, 8)` and places it below the 5 rows of `stats`: the result has shape `(6, 8)`. ::: --- ## Splitting Arrays Splitting arrays is the process of dividing an array into multiple sub-arrays. This can be useful for various data manipulation tasks. Below, we'll explore some common methods for splitting arrays using NumPy. ### Split The `np.split()` function splits an array into multiple sub-arrays along a specified axis. ```{python} # Create an array arr = np.array([1, 2, 3, 4, 5, 6]) # Split the array into 3 equal parts split_array = np.split(arr, 3) for idx in range(len(split_array)): print(f"Split {idx + 1}:\n{split_array[idx]}") ``` ### Array Split The `np.array_split()` function allows you to split an array into unequal parts if needed. For an array of length `l` that should be split into `n` sections, it returns `l % n` sub-arrays of size `l//n + 1` and the rest of size `l//n`. ```{python} # Split the array into 4 parts (unequal) split_array_unequal = np.array_split(arr, 4) for idx in range(len(split_array_unequal)): print(f"Split {idx + 1}:\n{split_array_unequal[idx]}") ``` ### Horizontal Split The `np.hsplit()` function splits an array horizontally (column-wise). ```{python} # Create a 2D array array_2d = np.array([[1, 2, 3], [4, 5, 6]]) # Split the 2D array into 3 columns hsplit_array = np.hsplit(array_2d, 3) for idx in range(len(hsplit_array)): print(f"Split {idx + 1}:\n{hsplit_array[idx]}") ``` ### Vertical Split The `np.vsplit()` function splits an array vertically (row-wise). ```{python} # Split the 2D array into 2 rows vsplit_array = np.vsplit(array_2d, 2) for idx in range(len(vsplit_array)): print(f"Split {idx + 1}:\n{vsplit_array[idx]}") ``` ### Splitting at Given Positions Instead of a number of equal parts, all the split functions also accept a **list of positions** where to cut. The three groups of statistics have different sizes: attack (goals and assists), defence (clean sheets and defensive contributions), and discipline (yellow cards, red cards, own goals, and penalties missed). We cut the columns before index 2 and before index 4: ```{python} attack, defence, discipline = np.hsplit(stats, [2, 4]) print(attack.shape, defence.shape, discipline.shape) ``` In FPL, a midfielder gets 5 points for each goal, 3 for each assist, 1 for each clean sheet, and 2 for each match with defensive contributions. He loses 1 point for each yellow card, 3 for each red card, 2 for each own goal, and 2 for each penalty missed. We store the points in a one-dimensional array, in the same order as the columns of `stats`: ```{python} # FPL points of a midfielder for each statistic points = np.array([5, 3, 1, 2, -1, -3, -2, -2]) ``` ::: {.callout-note} ## Exercise: Split the Points Split `points` into three arrays, one for each group of statistics: attack (the first two values), defence (the next two), and discipline (the last four). Then check your split by computing the sum of each group: one goal and one assist are worth 8 points, one clean sheet and one defensive contribution are worth 3 points, and one of each disciplinary event costs 8 points. ::: ::: {.callout-tip collapse="true"} ## Solution ```{python} attack_p, defence_p, discipline_p = np.split(points, [2, 4]) print(np.sum(attack_p), np.sum(defence_p), np.sum(discipline_p)) ``` `points` is one-dimensional, so we use `np.split()` (or `np.hsplit()`, which behaves the same for 1D arrays). ::: --- ## Challenge: One Array for the Whole League A friend in your mini-league sent you the statistics of the 5 midfielders in their squad **stat by stat**, as one long array: the 5 goals, then the 5 assists, and so on. They follow the same order of statistics as our team. ```{python} rival_export = np.array([ 11, 3, 7, 1, 5, # goals 6, 5, 9, 2, 8, # assists 10, 12, 7, 13, 9, # clean sheets 3, 16, 2, 25, 6, # defensive contributions 2, 5, 1, 8, 3, # yellow cards 0, 0, 0, 1, 0, # red cards 0, 1, 0, 0, 0, # own goals 1, 0, 0, 0, 1, # penalties missed ]) ``` Without writing any loop: 1. Turn `rival_export` into an array with one row per player and one column per statistic. 2. Join it with our `stats` into a single array, `league`, with our players first. Print its shape. 3. Split `league` into the three groups of statistics, and compute the fantasy points of each player in each group (use `points`). 4. Compute the total fantasy points of every player. 5. Build a report with two columns: the player number (1 to 10) and the total points. ::: {.callout-tip collapse="true"} ## Possible solution ```{python} # 1. One row per statistic first, then transpose rival = rival_export.reshape(8, 5).T print(rival) # 2. Join along the rows league = np.concatenate((stats, rival), axis=0) # or np.vstack() print(league.shape) # 3. Split the columns at the start of each group attack, defence, discipline = np.hsplit(league, [2, 4]) attack_p, defence_p, discipline_p = np.split(points, [2, 4]) print(np.sum(attack * attack_p, axis=1)) print(np.sum(defence * defence_p, axis=1)) print(np.sum(discipline * discipline_p, axis=1)) # 4. Total fantasy points total = league @ points print(total) # 5. Report report = np.column_stack((np.arange(1, 11), total)) print(report) ``` Note that `rival_export.reshape(5, 8)` would also give the right shape, but with the statistics mixed up. ::: --- ## Exam-Style Practice The midterm and the final exam have the same kinds of exercises you have been solving: a debugging exercise and an open problem. Try these two without AI, as you would in the exam, and open the solution only when you are done. Both come from the course of Manuel Morante. ::: {.callout-note} ## Debugging Exercise: Padding a Matrix The following function should add a border of a given number around a matrix, like the function `np.pad()`. It adds one row above, one row below, one column on the left, and one column on the right. This can be useful when every element needs neighbors in all directions. However, it raises this error: ```text ValueError: all the input arrays must have same number of dimensions, but the array at index 0 has 1 dimension(s) and the array at index 1 has 2 dimension(s) ``` Explain the error and fix the function. Does your fix work for a matrix that is not square? ```python def manual_pad(A, number): # Create rows and columns to add around the current array rows = np.ones(A.shape[0]) * number columns = np.ones(A.shape[1]) * number # Stack vertically and horizontally the new rows and columns C = np.vstack((rows, A, rows)) D = np.hstack((columns, C, columns)) return D A = np.array([[1, 0, -1], [-2, 1, 2], [0, -3, 1]]) print(manual_pad(A, 1)) ``` ::: ::: {.callout-tip collapse="true"} ## Solution There are three problems, and each one is about shapes: 1. `np.vstack()` accepts the 1D `rows` because it treats a vector as a row, so the vertical stack works. 2. `np.hstack()` needs arrays with the same number of dimensions. `columns` is 1D, and `C` is 2D. This is the error. A column must have the shape `(n, 1)`, not `(n,)`. 3. The new columns must be as long as the array **after** adding the rows, which is `C.shape[0]`, not `A.shape[1]`. Also, a new row must have `A.shape[1]` values, not `A.shape[0]`. The code only works by chance when the matrix is square, and a rectangular matrix would fail with a different error. ```{python} def manual_pad(A, number): # One row with as many values as columns, shape (1, columns) rows = np.full((1, A.shape[1]), number) C = np.vstack((rows, A, rows)) # One column with as many values as the new number of rows, shape (rows, 1) columns = np.full((C.shape[0], 1), number) return np.hstack((columns, C, columns)) A = np.array([[1, 0, -1], [-2, 1, 2], [0, -3, 1]]) R = np.array([[1, 2, 3, 4], [5, 6, 7, 8]]) # not square print(manual_pad(A, 1)) print(manual_pad(R, 0)) ``` NumPy already has this function. Check that your version agrees with it: ```{python} print(np.array_equal(manual_pad(R, 7), np.pad(R, 1, constant_values=7))) ``` Knowing that `np.pad()` exists is what matters in the exam, because you will find it in the cheat sheet. Understanding the shapes is what lets you use it, or write it yourself when it is not there. ::: ::: {.callout-note} ## Open Problem: Vectorized Evaluation of Mathematical Functions Your mathematics teacher asks you to apply several functions to a range of numbers. The challenge is to use a **single array** for the whole range of input values. Do not create a separate input array for each function. 1. Create a 1D array with **21 equally spaced values from 0.5 to 4**, including both endpoints. Its length is divisible by 3. 2. Split it into three equal sub-arrays. Evaluate the following functions on the corresponding sub-array, with vectorized operations and no loops: - $f_1(x) = x^3 - 3x + 2$ on the first sub-array. - $f_2(x) = \ln(x) + \sqrt{x}$ on the second sub-array. - $f_3(x) = e^{x/2}$ on the third sub-array. 3. For each function, find the positions where its result is greater than `2`. A position is the index inside its own sub-array, from 0 to 6. 4. Build a 2D array with three columns, one for each function, in this order. Each column has the positions of its function. The columns may have different lengths, so fill the shorter ones with `np.nan`. ::: ::: {.callout-tip collapse="true"} ## Solution ```{python} x = np.linspace(0.5, 4, 21) x1, x2, x3 = np.split(x, 3) f1 = x1**3 - 3*x1 + 2 f2 = np.log(x2) + np.sqrt(x2) f3 = np.exp(x3 / 2) # One array of positions for each function positions = [np.where(f > 2)[0] for f in (f1, f2, f3)] print(positions) ``` The columns have different lengths, so we first create the final array full of `np.nan` and then fill each column with slicing. The number of rows is the length of the longest column: ```{python} n_rows = max(len(p) for p in positions) result = np.full((n_rows, 3), np.nan) for column, p in enumerate(positions): result[:len(p), column] = p print(result) ``` On this range, $f_1$ never exceeds `2` (its largest value is about `1.07`), so `np.where` returns an empty array and the first column is all `np.nan`. A correct solution must not fail in that case. The loop above only assembles the result: the functions themselves are applied without loops. ::: --- ## Exam-Style Practice The midterm and the final exam have the same kinds of exercises you have been solving: a debugging exercise and an open problem. Try each one on your own, without AI, as you would in the exam, before opening its collapsed solution. ::: {.callout-note} ## Debugging: Padding a Matrix *Adapted from the exam-style exercises of Manuel Morante.* The following function should add a border of a given number around a matrix, like the function `np.pad()`. It adds one row above, one row below, one column on the left, and one column on the right. This is useful when every element needs neighbors in all directions. However, it raises this error: ```text ValueError: all the input arrays must have same number of dimensions, but the array at index 0 has 1 dimension(s) and the array at index 1 has 2 dimension(s) ``` Explain the error and fix the function. Does your fix work for a matrix that is not square? ```{python} #| eval: false def manual_pad(A, number): # Create rows and columns to add around the current array rows = np.ones(A.shape[0]) * number columns = np.ones(A.shape[1]) * number # Stack vertically and horizontally the new rows and columns C = np.vstack((rows, A, rows)) D = np.hstack((columns, C, columns)) return D A = np.array([[1, 0, -1], [-2, 1, 2], [0, -3, 1]]) print(manual_pad(A, 1)) ``` ::: ::: {.callout-tip collapse="true"} ## Solution There are three problems, and each one is about shapes: 1. `np.vstack()` accepts the 1D `rows` because it treats a vector as a row, so the vertical stack works. 2. `np.hstack()` needs arrays with the same number of dimensions. `columns` is 1D, and `C` is 2D. This is the error. A column must have the shape `(n, 1)`, not `(n,)`. 3. The sizes are wrong as well. A new row must have `A.shape[1]` values (one per column), not `A.shape[0]`, and a new column must be as long as the array **after** adding the rows, which is `C.shape[0]`. The original code only has matching sizes for a square matrix, so a rectangular one would still fail after fixing the dimensions. ```{python} def manual_pad(A, number): # One row with as many values as columns, shape (1, columns) rows = np.full((1, A.shape[1]), number) C = np.vstack((rows, A, rows)) # One column with as many values as the new number of rows, shape (rows, 1) columns = np.full((C.shape[0], 1), number) return np.hstack((columns, C, columns)) A = np.array([[1, 0, -1], [-2, 1, 2], [0, -3, 1]]) R = np.array([[1, 2, 3, 4], [5, 6, 7, 8]]) # not square print(manual_pad(A, 1)) print(manual_pad(R, 0)) ``` NumPy already has this function. Check that your version agrees with it: ```{python} print(np.array_equal(manual_pad(R, 7), np.pad(R, 1, constant_values=7))) ``` Knowing that `np.pad()` exists is what matters in the exam, because you will find it in the cheat sheet. Understanding the shapes is what lets you use it, or write it yourself when it is not there. ::: ::: {.callout-note} ## Open Problem: Vectorized Evaluation of Mathematical Functions *Adapted from the exam-style exercises of Manuel Morante.* Your mathematics teacher asks you to apply several functions to a range of numbers. The challenge is to use a **single array** for the whole range of input values. Do not create a separate input array for each function. 1. Create a 1D array with **21 equally spaced values from 0.5 to 4**, including both endpoints. Its length is divisible by 3. 2. Split it into three equal sub-arrays. Evaluate the following functions on the corresponding sub-array, with vectorized operations and no loops to apply them: - $f_1(x) = x^3 - 3x + 2$ on the first sub-array. - $f_2(x) = \ln(x) + \sqrt{x}$ on the second sub-array. - $f_3(x) = e^{x/2}$ on the third sub-array. 3. For each function, find the positions where its result is greater than `2`. A position is the index inside its own sub-array, from 0 to 6. 4. Build a 2D array with three columns, one for each function, in this order. Each column has the positions of its function. The columns may have different lengths, so fill the shorter ones with `np.nan`. *Expected result:* ```text [[nan 1. 0.] [nan 2. 1.] [nan 3. 2.] [nan 4. 3.] [nan 5. 4.] [nan 6. 5.] [nan nan 6.]] ``` ::: ::: {.callout-tip collapse="true"} ## Solution ```{python} x = np.linspace(0.5, 4, 21) x1, x2, x3 = np.split(x, 3) f1 = x1**3 - 3*x1 + 2 f2 = np.log(x2) + np.sqrt(x2) f3 = np.exp(x3 / 2) # One array of positions for each function positions = [np.where(f > 2)[0] for f in (f1, f2, f3)] print(positions) ``` The columns have different lengths, so we first create the final array full of `np.nan` and then fill each column with slicing. The number of rows is the length of the longest column: ```{python} n_rows = max(len(p) for p in positions) result = np.full((n_rows, 3), np.nan) for column, p in enumerate(positions): result[:len(p), column] = p print(result) ``` On this range, $f_1$ never exceeds `2` (its largest value is about `1.07`), so `np.where` returns an empty array and the first column is all `np.nan`. A correct solution must not fail in that case. The loop above only assembles the result: the functions themselves are applied without loops. ::: --- ## Summary You do not need to memorize every joining and splitting function. Focus on the shape of each input, the shape you want as output, and the axis along which the operation acts. - `reshape()` changes the arrangement of values but not their count. - `flatten()` returns a one-dimensional copy. - `concatenate()` joins arrays along an existing axis. - `stack()`, `vstack()`, and `hstack()` provide common arrangements. - `split()` requires equal pieces; `array_split()` also permits unequal pieces. Before running any joining operation, predict the output shape. Then explain why all non-joining dimensions must match. --- # Homework and Review Work through the homework exercises available here. [Homework](homework_numpy.qmd#array-manipulation){.btn .btn-primary} ---