All content here is under a Creative Commons Attribution CC-BY 4.0 and all source code is released under a BSD-2 clause license. Parts of these materials were inspired by https://github.com/engineersCode/EngComp/ (CC-BY 4.0), L.A. Barba, N.C. Clementi.
Please reuse, remix, revise, and reshare this content in any way, keeping this notice.
In the prior module you created vectors, matrices and arrays, with numbers in them. Now you get to perform calculations.
You will learn these by discovery: copy and paste the code, but look at the code first. Often there are other interesting ideas from prior modules, or even other new features we will learn about later.
Start a new version controlled repository for your work. Put each section in a new Python script file. Commit your work at the end of each section.
Learn by discovering...
import numpy as np
# Create an array with 3 rows and 4 columns
x = np.array([[1, 2, 3, 4],
[5, 6, 7, 8],
[9, 10, 11, 12]])
# Access individual values:
print(x[0, 0]) # first
print(x[2, 3]) # last
# Get vector "slices"
print(x[0, :]) # Row 1
print(x[1, :]) # Row 2
print(...) # Column 3
# Get sub-arrays:
print(x[1:3, 1:3]) # what section is this?
print(x[1:, 1:3]) # and this?
print(x[1:, 1:]) # and this?
# Negative indexing is also allowed, like we saw with lists.
print(x[-2:, :]) # The last 2 rows
You can also overwrite entries in the array:
# Create an array with 3 rows and 4 columns
x = np.array([[1.1, 2, 3, 4],
[5, 6, 7, 8],
[9, 10, 11, 12]])
# Change individual values:
x[0, 0] = 0
x[1, 2] = 3
print(x)
# Change an entire slice in row 2:
x[1, :] = 2
# Change the last column to be all -4
x[:, -1] = -4
# Reset the entire matrix to infinty!
x[:, :] = np.inf
Once we have created an array - see the prior notebook - we are then ready to actually use them for calculations!
Let us consider these calculations:
Commit your prior work, and start a new Python script for this section.
Create 2 matrices, of the same shape, and let's get started:
import numpy as np
A = np.ones(shape=(5,5))
B = np.ones(shape=(5,5))
print('A = \n{}\n\nB = \n{}'.format(A, B))
print(A + B)
The +
operation on two arrays is actually just a convenience. The actual function in NumPy which is being called to do the work is the np.add(...)
function.
np.add(A, B) # does exactly the same
Similarly, we have the -
and np.subtract()
functions that serve the same purpose:
print(A - B)
print(np.subtract(A, B))
print(np.add(A, -B))
These are called *element-by-element operations*. That means, NumPy performs the operation of addition on each corresponding element in the arrays A
and B
and then repeats that entry-by-entry. This is also called elementwise in NumPy's documentation.
NumPy will also allow you take shortcuts. Imagine that you want to subtract the value of 3.5 from every entry in matrix A
. You do not first need to create a matrix with the same shape as A
with all entries containing the values of 3.5, and then go to subract that.
There is a shortcut:
# Also does element-by-element calculations
print(A - 3.5)
Multiplication and division can also be done element-by-element using the *
and /
operators.
Try these to learn something new:
Create a matrix with the numbers 1 to 25 inclusive,
reshape
d into 5 rows and 5 columns.
- Divide each entry by 10.
- Now multiply each entry by
np.pi
$\approx 3.14159$- Start from the same matrix, and divide each entry by
np.e
$\approx 2.71828$- Try dividing by zero, in regular Python (i.e. not with NumPy):
5.6 / 0
- Now try dividing a NumPy matrix by zero: why do you get a different result? Does the NumPy result make more sense, or less sense, than regular division by zero in Python?
*Important note for MATLAB users*
Matlab is unique among programming languages for overloading the *
operator for matrix multiplication. In MATLAB that requires the matrices to the left and right of *
to have a consistent shape. We will see the NumPy equivalent of regular matrix multiplication further down.
There are other elementwise operations that can be done on matrices. These often involve raising the individual matrix elements to a certain power, or taking a square root (which is the same as raising a number to the power of $0.5$), or calculating the logarithm.
Let's try it out interactively.
- Use the
**
operation to raise a matrix to a power- Use the
.square()
function- Use the
.power()
function- Use the
.sqrt()
function- Verify that
**(0.5)
gives the same values as the.sqrt()
function
Extend your knowledge: try the above on a vector or matrix that contains values which are negative, zero and positive. What do you notice about the square root of a negative number?
A wide variety of mathematical functions are possible to be calculated on NumPy arrays. See the full list in the NumPy documentation.
You will self-discover these function by running the code below.
Try the following on this matrix:
# A 4x4 matrix with positive and negative values:
values = np.reshape(np.linspace(-2, +2, 16), (4, 4))
print(values)
- Try using the standard trigonometric functions
np.sin(...)
,np.tan(...)
, etc on that matrix.- Rounding each value to the closest integer. Use the
np.around()
function.Which way do negative values round?
- Round each value to a certain number of
decimals
; try rounding to 1 decimal place.Are the results what you expect?
- Similar to rounding: try the
np.floor(...)
andnp.ceil(...)
: what is the difference between the floor and the ceiling? Hint: read the documentation forfloor
andceil
.Watch for negative values: what do you expect
np.floor(3.5)
to do? Andnp.floor(-3.5)
- Logarithms and exponents are also part of the standard calculations we expect to do with matrices using the
np.log(...)
andnp.exp(...)
functions. Recall that $\log(\exp(x)) = x$.But this might be surprising to you:
exponent = np.exp(values)
print(exponent)
recovered = np.log(exponent)
print(recovered)
print('-----')
print(recovered - values)
> Does ``recovered`` match the original ``values`` matrix? It should: we first took the exponent, then the logarithm. This subtraction should be a matrix of all zeros.
>
> The last matrix in your printout above should be all zeros, but is not exactly equal to zero (it is very, very close to zero though).
>
> To test that we can use the ``np.isclose(...)`` function. It is another elementwise function that you can add to your toolbox. It tests if the entries in an array are close to another:
```python
np.isclose(recovered - values, 0)
This problem uses all the skills you have learned so far: creating NumPy arrays from a list, and doing calculations on them.
Back in module 1 we saw that $$n! \approx \sqrt{2 \pi n} \, \cdot \, n^n \, \cdot \, e^{-n} $$
is what is known as the Stirling approximation for the factorial of an integer value $n$.
N = 15
# 1. Create a vector of values first: 1, 2, ... N
n_values = np.linspace(...)
# 2. Find Stirling's approximation for the factorial
approximate = np.sqrt(...) * ....
# 3. True value of n! (to compare with later)
import math
true_values = []
for n in np.arange(N):
# Use math.factorial(n) to calculate true value
# Append that to the list
...
# 4. Convert list to a NumPy array
# 5. Finally, calculate and show the percentage relative error
error = ...
percentage_error = ...
*Observe:* how the problem above (to discover the approximation error) was broken down into 5 sub-steps. Try doing that for your own programs. Create a new file and break the problem down:
then turn those lines into comments, and get started coding.
<img src="images/general/Crystal_Clear_action_db_commit.png" style="width: 100px ; float:left"/> <br><br> ***Have you committed your script to your repo? Also push it to the remote server.***
You have just obtained all the data in your matrix, and now you wish to find the largest, or smallest value in the entire matrix. Or perhaps the average of all numbers. These are global operations.
import numpy as np
rnd = np.array([[ 7, 3, 11, 12, 2],
[10, 13, 8, 8, 2],
[ 3, 13, 6, 2, 3],
[ 5, 3, 9, 2, 6]])
print('The matrix is:\n{}'.format(rnd))
max_value = np.amax(rnd)
print('The maximum value is {}'.format(max_value))
min_value = np.amin(rnd)
print('The minimum value is {}'.format(min_value))
# Note here how strings can be split across lines:
print(('The mean value is {}, while the median '
'is {}').format(np.mean(rnd), np.median(rnd)))
The np.amax(...)
and np.amin(...)
functions will find the maximum or minimum in the entire array: looking at every single entry in all dimensions.
The NumPy library will internally unfold, or flatten the array into a single long vector. Take a look at what that looks like when you use the .flatten(...)
method on the array: rnd.flatten()
. It works by going from column to column, down each row:
print(rnd.flatten())
[ 7 3 11 12 2 10 13 8 8 2 3 13 6 2 3 5 3 9 2 6]
This is actually how the data is stored internally in the computer's memory.
The reason we point this array.flatten(...)
function out is because sometimes knowing what the maximum value is is only half the work. The other half is knowing where that maximum value is. For that we have the np.argmax(...)
function.
Try this code:
max_position = np.argmax(rnd)
print(('The maximum value is in position {} of the '
'flattened array').format(max_position))
Verify that that is actually the case, using the space below to run the code. The max_position
is the position in the flatten
ed array, which is why you should know about the array.flatten(...)
function.
Another tip: notice that the maximum value of 13
appears twice. To which one does max_position
refer?
Other operations on the entire matrix, which reduce the entire matrix down to a single number. Try these out:
# Related to maximum and minimim is `ptp`. Use
# the help function to see what this does.
rnd.ptp()
rnd.mean()
rnd.prod() # what does this do?
rnd.std()
rnd.sum()
rnd.trace() # sum of the elements on the diagonal
rnd.var()
Above you learned about elementwise operations. In other words, NumPy performed the mathematical calculation on every element (entry) in the array. You also saw functions which summarize the matrix by 1 single value (like the mean or standard deviation).
Sometimes we need calculations that work on every row, or column, of an array. We will cover:
We will talk about matrices in these examples, but the operations can also be applied to multi-dimensional arrays, with 3, 4, or more dimensions.
We also introduce 2 important terms: *
axis
* and *inplace
*, both of which you will regularly see in the NumPy documentation, but also later in Pandas.
Before you get started, commit prior work, and start this section in a new script.
Think of a matrix containing the daily temperatures per city; one city per column, and every row is a day of the year.
City 1 | City 2 | City 3 | City 4 | |
---|---|---|---|---|
Day 1 | 7 | 9 | 12 | 10 |
Day 2 | 1 | 4 | 5 | 2 |
Day 3 | -3 | 1 | -2 | -3 |
Day 4 | -2 | -1 | -2 | -2 |
Day 5 | -3 | -1 | -2 | -4 |
For this we also use the np.amax(array, axis=...)
or np.amin(array, axis=...)
function, which you saw above. But, now you must specify, as a second input, along which *axis
* you want that extreme value to be calculated.
Try out the code below. Note that the array is created from a list-of-lists.
import numpy as np
temps = np.array([[7, 9, 12, 10],
[1, 4, 5, 2],
[-3, 1, -2, -3],
[-2, -1, -2, -2],
[-3, -1, -2, -4]])
print(('The temperatures are given one column per '
'city, each row is a daily average '
'temperature:\n{}').format(temps))
max_value_0 = np.amax(temps, axis=0)
print(('The maximum value along axis 0 (row-wise, '
'per city for all days) is {}').format(max_value_0))
max_value_1 = np.amax(temps, axis=1)
print(('The maximum value along axis 1 (column-wise, '
'per day for all cities) is {}').format(max_value_1))
# Notice the above output is 'flatten'ed and returned as a row,
# instead of a column, as you might hope for. We can use
# the `keepdims` input though:
max_value_1_col = np.amax(temps, axis=1, keepdims=True)
print(('The maximum value along axis 1 (column-wise, '
'per day for all cities) is\n{}').format(max_value_1_col))
You can visually verify that the maximum values returned are what you expected.
Now try it for the minimum values of the matrix:
Many functions in NumPy take *axis
* as in input argument, including the np.mean(...)
, np.std(...)
and other functions you saw above.
Complete this code:
average_temperature_per_city = ...
print(('The average temperature for each city '
'is: {}').format(average_temperature_per_city))
average_temperature_per_day = ...
print(('The average temperature per day for all '
'cities: {}').format(average_temperature_per_day))
std_per_city = ...
print(('The standard deviation of the temperature '
'is {}').format(std_per_city))
Just like finding the minimum or maximum value in every column, you might also be interested in sorting each column.
We want to sort every column independently of the others. In other words every column will be sorted from low to high by the end. This is in contrast to sorting based on one column, and the other rows follow with. We will see that coming later.
We have seen the *axis
* input several times now, and here we will use it again to indicate which axis we would like to sort along.
import numpy as np
temps = np.array([[7, 9, 12, 10],
[1, 4, 5, 2],
[-3, 1, -2, -3],
[-2, -1, -2, -2],
[-3, -1, -2, -4]])
sorted_columns = np.sort(temps, axis=0)
print('The temperatures, sorted in each column, along axis 0: \n{}'.format(sorted_columns))
sorted_rows = np.sort(temps, axis=1)
print('The temperatures, sorted in each row, along axis 1: \n{}'.format(sorted_rows))
print('Verify the original matrix is untouched:\n{}'.format(temps))
The sort takes place and the result is set to a new variable. This is not efficient, especially if the input_array to be sorted is really large. It means that a copy of the data is made, using up memory, and computer time; and only then does the sort take place on the copy of the data.
It is possible to simply sort the original array. This is called *in-place sorting. You will see that terminology in several places in NumPy's documentation: in-place*. It means that the original matrix is used, calculated on, and the result is left in the same variable as the original matrix.
You perform an in-place sort as follows:
input_array.sort(axis=...)
compared with a sort that is assigned to a new variable:
output_array = np.sort(input_array, axis=...)
By definition an in-place operation means there is no need to assign the result to another variable as output. The second variation with np.sort(input_array, ...)
follows good practice, where a function does not modify its inputs. So therefore the output_array
is required.
Let's see what the implication of in-place sorting is:
import numpy as np
temps = np.array([[7, 9, 12, 10],
[1, 4, 5, 2],
[-3, 1, -2, -3],
[-2, -1, -2, -2],
[-3, -1, -2, -4]])
# In-place sort. We don't need to use `output=` but let's see what happens.
# The in-place sort result is in the original variable "temps"
output = temps.sort(axis=0)
print('The sorted values along axis 0: \n{}'.format(temps))
print('Out of curiosity, the value of "output" is: {}'.format(output))
# So you can simply say:
temps.sort(axis=0)
# and the result will be sorted in the original variable.
Note the above behaviour is consistent with what we saw before in regular Python.
If the values in the rows or columns represent some property, such as weight in a container, then it might be interesting to calculate the cumulative value.
We will create a matrix where the values in column 0 will all be ones. Values in column 1, 2, 3 and 4 will represent the weight (kilograms), added to 4 containers. Each row represents 1 minute. Column 0 is a counter.
The goal: Find the point in time when the weight in the container just exceeds 100kg. You will see why we have a counter as column 1.
import numpy as np
n = 20
k = 4
counter = np.ones(shape=(n, 1))
weights = np.random.randint(low=4, high=9, size=(n, 4))
# Stack the the objects side-by-side (h=horizontal)
weight_matrix = np.hstack((counter, weights))
print('The weight added to the {} containers every minute [ignore first column]:\n{}'.format(k, weight_matrix))
accumulation = np.cumsum(weight_matrix, axis=0)
print('The cumulative weight over time is:\n{}'.format(accumulation))
The methods to calculate something globally, or something on just a row or a column are the same methods: just add the axis
keyword:
np.sum
is over the entire arraynp.sum(axis=...)
is over the specified axis.Operation | Description |
---|---|
np.sum |
Sum of entries |
np.prod |
Product (multiplication) of all entries |
np.mean |
Arithmetic average |
np.std |
Standard deviation |
np.var |
Variance |
np.min |
Minimum value |
np.max |
Maximum value |
np.argmin |
Index of the minimum |
np.argmax |
Index of the maximum |
np.median |
Median value |
Use the skills just learned above regarding calculations, and in the prior module regarding random numbers, to solve the following 2 problems.
Complete this challenge in a new script. Have you pushed your prior work?
Simulate $n=[40, 400, 4000, 40000, 4\times 10^6]$ coin tosses in Python, in a single NumPy vector. For each one of these 5 scenarios, calculate the average of the vector.
The first column is the station number, the next is the date, and the third column is the temperature, measured in units of °C.
filename = 'KNMI-Homogenized-temperature-monthly-average-1901-2019.txt' full_filename = os.path.join(...) weather = np.loadtxt(full_filename, skiprows=27)
* What is the lowest temperature recorded over that time period? And the highest?
* What is the standard deviation of the temperature values?
* What is the average of the first 700 rows in the data set? And the last 700 rows?
* Repeat for the standard deviation as well.
* What is the average temperature since the year you were born? And prior?
In future modules we will plot and examine other statistics from this data file.
Matrices can be considered like spreadsheets containing rows and observations. But the real power of matrices comes from some of the calculations we can do with them.
We can multiply a matrix, called $\mathbf{A}$ that has 4 rows and 2 columns, with a second matrix, $\mathbf{B}$, which contains 2 rows and 3 columns. The product of $\mathbf{A}$ and $\mathbf{B}$, or also written as $\mathbf{AB}$ is another matrix, we have called it $\mathbf{X} = \mathbf{AB}$. That matrix $\mathbf{X}$ has 4 rows and 3 columns.
The value in first row and first column of $\mathbf{A}$ is given by the variable $a_{1,1}$; the value in row 1 and column 2 is $a_{1,2}$, and so on. The entry in row $i$ and column $j$ is therefore given by $a_{i,j}$. We use a similar notation for matrix $\mathbf{B}$ and $\mathbf{X}$. You can see these values in the figure shown here. Figure credit: Wikipedia.
$$\begin{eqnarray*}\mathbf{A}\mathbf{B} &=& \mathbf{X}\\ \begin{bmatrix} {a_{1,1}} & {a_{1,2}} \\ {a_{2,1}} & {a_{2,2}} \\ {a_{3,1}} & {a_{3,2}} \\ {a_{4,1}} & {a_{4,2}} \\ \end{bmatrix} \begin{bmatrix} {b_{1,1}} & {b_{1,2}} & {b_{1,3}} \\ {b_{2,1}} & {b_{2,2}} & {b_{2,3}} \\ \end{bmatrix} &=& \begin{bmatrix} x_{1,1} & x_{1,2} & x_{1,3} \\ x_{2,1} & x_{2,2} & x_{2,3} \\ x_{3,1} & x_{3,2} & x_{3,3} \\ x_{4,1} & x_{4,2} & x_{4,3} \\ \end{bmatrix}\end{eqnarray*} $$The values at the intersection of $\mathbf{X}$ are marked with circles:
$\begin{align*} x_{1,2} & = {a_{1,1}}{b_{1,2}} + {a_{1,2}}{b_{2,2}} \\ x_{3,3} & = {a_{3,1}}{b_{1,3}} + {a_{3,2}}{b_{2,3}}\\ \end{align*}$
Notice how the dimensions are consistent: a vector of 1 row and 2 columns $[{a_{1,1}} \,\,\, {a_{1,2}}]$ is multiplied element-by-element with $[{b_{1,2}}\,\,\, {b_{2,2}}]$, and then the multiplied values are added. That process is repeated over and over, until all the values in the new matrix $\mathbf{X}$ is formed.
There is also a pattern: the value that goes into row 1 and column 2 of $\mathbf{X}$ comes from multiplying row 1 of $\mathbf{A}$ and column 2 of $\mathbf{B}$. This form of multiplying the vectors element-by-element and then adding up the result is called the *dot product*. This explains why the function to do matrix multiplication in NumPy is called np.dot(...)
.
Just a quick mention: you cannot just multiply any two matrices $\mathbf{A}$ and $\mathbf{B}$. We do require that:
There is a lot of multiplying and adding happening here; 12 dot products. This is certainly not something you want to do by hand regularly. So ... let's try it with NumPy.
import numpy as np
A = np.array([[8, 7], [6, 5], [4, 3], [2, 1]])
B = np.array([[1, 2, 3], [4, 5, 6]])
print('The shape of matrix A is: {}\n'.format(A.shape))
print('The shape of matrix B is: {}\n'.format(B.shape))
X = np.dot(A, B)
print('The matrix multiplication of AB = X \n{}'.format(X))
# The value of X can also be calculated from:
X = A.dot(B)
print(('The value in position (2,3) of X comes from row 2 '
'of A [{}] and column 3 of B [{}] which corresponds '
'to 6*3 + 5*6 = {}').format(A[1,:], B[:,2], 6*3 + 5*6))
@
operator for matrix multiplication¶Recently (in 2014) Python introduced a new operator: @
specifically defined for matrix multiplication. Depending on if your Python version is up to date, you do not need to write np.dot(A, B)
anymore, but rather, more compactly: A @ B
for matrix multiplication. You will still see np.dot(A, B)
in older NumPy code, but A @ B
is much cleaner.
# Requires a recent Python version
import numpy as np
A = np.array([[8, 7], [6, 5], [4, 3], [2, 1]])
B = np.array([[1, 2, 3], [4, 5, 6]])
X = A @ B
print('Compact matrix multiplication: A@B = \n{}'.format(X))
# which is the same result as given above.
*Try this important exercise:*
Let $\mathbf{A}$ be a matrix with 4 rows and 3 columns. Let $\mathbf{B}$ be a matrix with 3 rows and 2 columns.
import numpy as np A = np.array([[8, 7, 3], [6, 5, 4], [4, 3, 5], [2, 1, 6]]) B = np.array([[1, 2], [4, 5], [7, 8]])
>
>Verify if $\mathbf{AB}$ is equal to $\mathbf{BA}$. *Hint*: we already know from the shapes of $\mathbf{A}$ and $\mathbf{B}$ that this is not true. Remember the rule about the number of rows and columns?
>
> Now try it when the number of rows in $\mathbf{B}$ do match the number of columns in $\mathbf{A}$:
> ```python
A = np.array([[8, 7], [6, 5]])
B = np.array([[1, 2], [3, 4]])
Is it true in general that $\mathbf{AB} = \mathbf{BA}$, or equivalently: $\mathbf{AB} - \mathbf{BA} = \mathbf{0}$?
In the above example you saw an error message when the dimensions of the matrix don't match. For example, you have two matrices:
and you want to multiply them to get $\mathbf{X}=\mathbf{AB}$ with with 4 rows and 2 columns.
Matrix $\mathbf{B}$ does is not in the right order: it has 2 rows and 3 columns, but it should have 3 rows and 2 columns
We first need to first *transpose* $\mathbf{B}$, which means to flip the rows and columns: rows become columns and columns become rows.
In this example once we have transposed matrix $\mathbf{B}$ it will have 3 rows now and 2 columns, and be ready to multiply. We therefore write: $\mathbf{X}=\mathbf{AB}^T$, or you sometimes also see $\mathbf{B}'$ to indicate the transpose.
There above you see the notation we use to indicate a transpose:
Transposing happens so frequently that there is shortcut for it in NumPy: B.T
will transpose matrix B
.
*Try it:*
import numpy as np
A = np.array([[8, 7, 3], [6, 5, 4], [4, 3, 5], [2, 1, 6]])
B = np.array([[1, 2, 3], [4, 5, 6]])
print('Matrix B has this shape: {}'.format(B.shape))
B_transpose = B.T
print('B transposed is:\n{}'.format(B_transpose))
# Stack the operations in one line!
print('B transposed has shape: {}'.format(B.T.shape))
We don't need to create an intermediate variable B_transpose
to hold the transposed value.
Finally, to calculate $\mathbf{AB}^T$ we can simply write: np.dot(A, B.T)
or even shorter A @ B.T
. Complete the code:
# (out loud you would say: "A times B transposed")
print('A times B transposed is {}'.format( ... ))
A square matrix (a matrix with the same number of rows and columns) can be *inverted*. If the matrix is $\mathbf{A}$, then the inverse is written as $\mathbf{A}^{-1}$. Inversion is related to the concept of division. We won't go into the mathematical details of how an inverse is calculated though.
The command for inverting is:
A = np.random.random((4, 4))
A_inv = np.linalg.inv(A)
Other linear algebra commands of interest in the linalg
section of NumPy are:
np.linalg.eig
for eigenvalues and eigenvectorsnp.linalg.eigvals
for the eigenvalues onlynp.linalg.svd
for the singular value decomposition (related to PCA)np.linalg.solve
for solving linear matrix equationsnp.linalg.solve
for solving linear matrix equationsType help(np.linalg.info)
for more details, and other linear algebra routines.
You can apply most concepts learned above by solving this problem.
If you build a least squares regression model with more than one variable you can assemble all you input data (predictors, columns used to make the prediction) in a single matrix $\mathbf{X}$. For example:
$$ y = b_0 + b_\text{P} x_\text{P} + b_\text{T} x_\text{T}$$has 2 prediction columns: $x_\text{P}$ and $x_\text{T}$; but there is also an intercept.
Assemble matrix $\mathbf{X}$ to have 3 columns: a column of ones, a column for $x_\text{P}$ and one for the values of $x_\text{T}$. The following data are from a designed experiment, where the 3 required columns are:
$$ \mathbf{X} = \begin{bmatrix} 1 & -1 & -1 \\ 1 & +1 & -1 \\ 1 & -1 & +1 \\ 1 & +1 & +1 \\ 1 & 0 & 0\\ 1 & 0 & 0 \end{bmatrix} $$*Hint*: There are many way to create the above matrix, but consider using the np.hstack
function to stack vectors next to each other (side-by-side).
To build a regression model you also need your vector of $\mathbf{y}$ values, the measured values: $$ \mathbf{y} = \begin{bmatrix} 55\\90 \\ 190 \\ 170 \\135 \\ 142 \end{bmatrix} $$
The regression coefficients, the 3 values of $\mathbf{b} = \left[ b_0, \, b_\text{P},\, b_\text{T}\right]$ can be found from the following linear algebra expression: $\mathbf{b} = (\mathbf{X}^T \mathbf{X})^{-1} \mathbf{X}^T \mathbf{y}$
- Calculate that vector $\mathbf{b}$ with 3 values. A good check that you got the correct values is if the first entry is equal to
np.average(y)
- Once you have the coefficients, you can make predictions from the regression model: $\hat{\mathbf{y}} = \mathbf{X} \mathbf{b}$. Check the dimensions of $\mathbf{X}$ and $\mathbf{b}$ to ensure they are correct.
- Lastly calculate the residuals, the prediction error = actual $\mathbf{y}$ minus predicted = $\mathbf{e} = \mathbf{y} - \hat{\mathbf{y}}$.
Wrap up this section by committing all your work. Have you used a good commit message? Push your work, to refer to later, but also as a backup.
*Feedback and comments about this worksheet?* Please provide any anonymous comments, feedback and tips.
# IGNORE this. Execute this cell to load the notebook's style sheet.
from IPython.core.display import HTML
css_file = './images/style.css'
HTML(open(css_file, "r").read())