session_1.Rmd

---
title: "Session 1"
author: "Antoine Vernet" 
subtitle: Statistics and Programming with R
output:
  ioslides_presentation:
    css: style.css
  beamer_presentation:
    slide_level: 2
---

```{r setup, echo = FALSE}
library(knitr)
library(rgl)
knit_hooks$set(webgl = hook_webgl)
```

## Lesson plan
Today we will cover:

- Fundamentals of linear algebra
    + Matrix algebra
- Fundamentals of calculus

You will have covered most of those things during highschool, so you should be pretty familiar with everything.

## Definitions

- Algebra: the study of mathematical symbol and the rules to manipulate them

- Linear algebra: the study of vector spaces and linear mapping between those spaces

- Calculus: the study of change

- Geometry: the study of shape, size and relative position

(Source: Wikipedia)


# Matrix algebra


## What is a matrix?

> - A matrix is a rectangular array of numbers.
> - An $m \times n$ matrix has $m$ rows and $n$ columns.

$${\bf A} = [a_{ij}] = 
\begin{bmatrix} 
    a_{11} & a_{12} & a_{13} & \dots  & a_{1n} \\
    a_{21} & a_{22} & a_{23} & \dots  & a_{2n} \\
    \vdots & \vdots & \vdots & \ddots & \vdots \\
    a_{m1} & a_{m2} & a_{m3} & \dots  & a_{mn}
\end{bmatrix}$$

> - $i$ is the index for rows, $j$ is the index for columns.

## Definition: square matrix

- A square matrix has the same number of rows and columns (i.e $n \times n$) 

$${\bf A} = [a_{ij}] = 
\begin{bmatrix} 
    a_{11} & a_{12} & \dots  & a_{1n} \\
    a_{21} & a_{22} & \dots  & a_{2n} \\
    \vdots & \vdots & \ddots & \vdots \\
    a_{n1} & a_{n2} & \dots  & a_{nn}
\end{bmatrix}$$

## Definition: row  and column vector

- A row vector is a $1 \times m$ matrix, it can be written as ${\bf x} \equiv (x_1, x_2, ..., x_m)$ 

- A column vector is a $n \times 1$ matrix, it can be written as $${\bf y} \equiv \begin{bmatrix}
        y_1 \\
        y_2 \\
        \vdots\\
        y_n \\
        \end{bmatrix}$$ 

## Definition: diagonal matrix

A square matrix ${\bf A}$ is a diagonal matrix if all the elements off diagonal are 0 ($a_{ij} = 0$, for all $i \neq j$).

$${\bf A} = [a_{ij}] = 
\begin{bmatrix} 
    a_{11} & 0 & \dots  & 0 \\
    0 & a_{22} & \dots  & 0 \\
    \vdots & \vdots & \ddots & \vdots \\
    0 & 0 & \dots  & a_{nn}
\end{bmatrix}$$


## Definition: identity and zero matrix

The $n \times n$ __identity matrix__, denoted ${\bf I}$ (or ${\bf I_n}$)  is the diagonal matrix with unity in each diagonal position and zero elsewhere:

$${\bf I} = 
\begin{bmatrix} 
    1 & 0 & \dots  & 0 \\
    0 & 1 & \dots  & 0 \\
    \vdots & \vdots & \ddots & \vdots \\
    0 & 0 & \dots  & 1
\end{bmatrix}$$

The $m \times n$ __zero matrix__, denoted ${\bf 0}$ is the $m \times n$ matrix with zero for all entries, this need not be a square matrix.

## Matrix addition

Two matrices of the same dimensions, for example two $m \times n$ matrices can be added element by element: ${\bf A} + {\bf B} = [a_{ij} + b_{ij}]$.

$${\bf A} + {\bf B} =
\begin{bmatrix} 
    a_{11} + b_{11} & a_{12} + b_{12} & a_{13} + b_{13} & \dots  & a_{1n} + b_{1n} \\
    a_{21} + b_{21} & a_{22} + b_{22} & a_{23} + b_{23} & \dots  & a_{2n} + b_{2n} \\
    \vdots & \vdots & \vdots & \ddots & \vdots \\
    a_{m1} + b_{m1} & a_{m2} + b_{m2} & a_{m3} + b_{m3} & \dots  & a_{mn} + b_{mn}
\end{bmatrix}$$

## Scalar multiplication

Given any real number $\gamma$, __scalar multiplication__ is defined as $\gamma{\bf A} \equiv [\gamma a_{ij}]$ or:

$$\gamma{\bf A} = [\gamma a_{ij}] = 
\begin{bmatrix} 
    \gamma a_{11} & \gamma a_{12} & \gamma a_{13} & \dots  & \gamma a_{1n} \\
    \gamma a_{21} & \gamma a_{22} & \gamma a_{23} & \dots  & \gamma a_{2n} \\
    \vdots & \vdots & \vdots & \ddots & \vdots \\
    \gamma a_{m1} & \gamma a_{m2} & \gamma a_{m3} & \dots  & \gamma a_{mn}
\end{bmatrix}$$

## Matrix multiplication

In order for matrix multiplication of matrix ${\bf A}$ and ${\bf B}$ to be possible, the column dimension of ${\bf A}$ must equal the row dimension of ${\bf B}$. If ${\bf A}$ is an $m \times n$ matrix and ${\bf B}$ an $n \times p$ matrix, matrix multiplication is defined as:
$$
{\bf AB} = \left[\sum\limits_{k = 1}^n a_{ik} b_{kj}\right]
$$

for example:

$$
\begin{bmatrix} 
    2 & -1 & 0 \\
    -4 & 1 & 0
\end{bmatrix}
\begin{bmatrix} 
    0 & 1 & 6 & 0 \\
    -1 & 2 & 0 & 1 \\
    3 & 0 & 0 & 0
\end{bmatrix} = \begin{bmatrix} 
    1 & 0 & 12 & -1 \\
    -1 & -2 & -24 & 1
\end{bmatrix}
$$


## Transpose

Let ${\bf A} = [a_{ij}]$ be an $m \times n$ matrix. The transpose of $\bf A$ denoted $\bf A'$ is the $n \times m$ matrix obtained by interchanging the rows and columns of $\bf A$, so ${\bf A'} \equiv [a_{ji}]$.
For example:
$$
{\bf A} = \begin{bmatrix} 
    1 & 0 & 12 & -1 \\
    -1 & -2 & -24 & 1
\end{bmatrix}, \, 
{\bf A'} = \begin{bmatrix} 
    1 & -1 \\
    0 & -2 \\
    12 & -24\\
    -1 & 1
\end{bmatrix}
$$

```{r, echo = FALSE, eval = FALSE} 
## Partitioned matrix multiplication
```

## Trace

The trace of a matrix is defined for square matrices. 
For any $n \times n$ matrix $\bf A$, the trace of this matrix, tr$\left(\bf A\right)$ is the sum of its diagonal elements: 

$$
\mathrm{tr}({\bf A}) = \sum\limits_{i = 1}^n a_{ii}
$$

## Inverse

An $n \times n$ matrix ${\bf A}$ has an inverse, denoted ${\bf A^{-1}}$, provided that ${\bf A^{-1}A} = {\bf I_n}$ and ${\bf AA^{-1}} = {\bf I_n}$. In this case ${\bf A}$ is said to be _invertible_ or _nonsingular_, otherwise it is said to be _noninvertible_ or _singular_.

## Linear independence

Let $\{{\bf x_1}, {\bf x_2}, ..., {\bf x_r}\}$ be a set of $n \times 1$ vectors. Those are linearly independent vectors if and only if, 
$${ \alpha_1 {\bf x_1} + \alpha_2 {\bf x_2} + ... + \alpha_r {\bf x_r} = {\bf 0} }$$
implies that ${\alpha_1 = \alpha_2 = ... = \alpha_r = 0}$

In case the set of vector is linearly dependent it means that at last one of the vectors can be written as a linear combination of the others.

## Rank

The rank of a matrix ${\bf A}$ is the maximum number of linearly independent columns of $\bf A$.

For example, 

$$
\begin{bmatrix} 
    1 & 20 \\
    5 & 100 \\
    3 & 60\\
    0 & 0
\end{bmatrix}
$$

> - can, at most, be of rank 2. Its rank is one because the second columns is 20 times the first one.

# Matrices in R

## A first example

```{r}
A <- matrix(data = c(1, 2, 3, 4, 5, 6), nrow = 2, byrow = FALSE)

A
```

Scalar multiplication

```{r}
A * 2
```

## A second example

```{r}
B <- matrix(data = c(1, 2, 3, 4, 5, 6), nrow = 3, byrow = FALSE)

B

A %*% B
```


# Systems of equations

## A simple system

$$
\begin{cases} x + y + z = 2  \quad \qquad\qquad(1)\\ 5x - 8y + 2z = 2 \ \ \quad \qquad(2)\\ -2x + 3y - 5z = -3 \qquad (3)\end{cases}
$$

Solve by addition (or method of elimination) and substitution.

First step, add $2 \times (1)$ to $(3)$ and add $-5 \times (1)$ to $(2)$.

$$
\begin{cases} x + y + z = 2 \\ 
              - 13y - 3z = -8\\
              5y - 3z = 1
\end{cases}
$$


## Solving a system of equations (cont'd)


Second step, it looks like we can eliminate $z$ easily by substracting $(3)$ in $(2)$

$$
\begin{cases} x + y + z = 2 \\ 
              y = \frac{-8 + 3 z}{-13}\\
              5y - 3z = 1
\end{cases}
$$

You can then replace in $(3)$

$$
5 \left(\frac{-8 + 3 z}{-13}\right) - 3z = 1
$$

$$
z = \frac{1}{2}
$$


## Solving a system of equations (cont'd 2)

$$
\begin{cases} x + y + z = 2 \\ 
              y = \frac{-8 + 3 z}{-13}\\
              z = \frac{1}{2}
\end{cases}
$$

You can then replace  $z$ in $(2)$

$$ y = \frac{1}{2}$$

And finally we can replace in $(1)$ to obtain

$$x = 1$$

## Solving a system of equation: the matrix method

It is very similar to the previous method

Write down the system of equation as a matrix

$$\begin{bmatrix}
1 & 1 & 1 & | & 2\\
5 & -8 & 2 & | & 2\\
-2 & 3 & -5 & | & -3
\end{bmatrix}
$$

What we want to do is to transform the previous matrix into the following one:

$$\begin{bmatrix}
1 & 0 & 0 & | & a\\
0 & 1 & 0 & | & b\\
0 & 0 & 1 & | & c
\end{bmatrix}
$$

where $a$, $b$ and $c$ are the roots of the system.

## Solving the system

The first row of the first column is already where we want it with $x = 1$. Let's make the second and third row of the first column 0 by substracting $5 \times$ row 1 into row 2 and adding $2 \times$ row 1 into row 3.

$$\begin{bmatrix}
1 & 1 & 1 & | & 2\\
0 & -13 & -3 & | & -8\\
0 & 5 & -3 & | & 1
\end{bmatrix}
$$

Let's now make the second row of the second column 1 by dividing row 2 by -13.

$$\begin{bmatrix}
1 & 1 & 1 & | & 2\\
0 & 1 & \frac{3}{13} & | & \frac{8}{13}\\
0 & 5 & -3 & | & 1
\end{bmatrix}
$$

## 

Now let's make the third row of the second column of row 1 and row 3 zeros using the new row 2 we created (we substract row 2 from row 1 and substract $5 \times$ row 2 from row 3).

$$\begin{bmatrix}
1 & 0 & \frac{10}{13} & | & \frac{18}{13}\\
0 & 1 & \frac{3}{13} & | & \frac{8}{13}\\
0 & 0 & -\frac{54}{13} & | & -\frac{27}{13}
\end{bmatrix}
$$

## 

Now let's make row 3 of column 3 1 by multiplying it by $-\frac{13}{54}$

$$\begin{bmatrix}
1 & 0 & \frac{10}{13} & | & \frac{18}{13}\\
0 & 1 & \frac{3}{13} & | & \frac{8}{13}\\
0 & 0 & 1 & | & \frac{1}{2}
\end{bmatrix}
$$

Now we can replace create a new row 1 and 2 eliminating the values in the third column by substracting $\frac{3}{13} \times$ row 3 from row 2 and substracting $\frac{10}{13} \times$ row 3 from row 1.

## Solution

$$\begin{bmatrix}
1 & 0 & 0 & | & 1\\
0 & 1 & 0 & | & \frac{1}{2}\\
0 & 0 & 1 & | & \frac{1}{2}
\end{bmatrix}
$$


## Solving equations in R


```{r}
a <- matrix(data = c(1, 1, 1, 5, -8, 2, -2, 3, -5), 
            nrow = 3, byrow = TRUE)

b <- c(2, 2, -3)

a

solve(a, b)

```

# Calculus

## Function

A function is a relation between a set of inputs and a set of permissible outputs, each input is related to exactly one output.

The notation is: 
$$f(x) = y
$$

where $f$ is the function, $x$ is the domain, $y$ is the image and $(x, y)$ is the graph of the function. Usually, we say that $f$ maps $x$ to $y$.


## Linear functions

A linear function is a map between 2 vector spaces. It preserves vector addition and scalar multiplication.
Most (all) of the functions we will be encountering in this class are linear functions.

They are very useful to study linear processes but also as approximations for non-linear processes (e.g. logistic regression).

```{r echo = FALSE}
## Vector addition

## Scalar multiplication
```

## Derivatives 

The first derivative of a univariate function $f: \mathbb{R} \mapsto \mathbb{R}$ is the function $f': \mathbb{R} \to \mathbb{R}$ defined as:
$$
f'(x) = \lim_{h\downarrow 0}\frac{f(x + h) - f(x)}{h}
$$

```{r, echo = FALSE, out.width = '50%', fig.retina = NULL, fig.align = 'center'}
knitr::include_graphics("./img/1024px-Lim-secant.png")
```


## Illustration {.flexbox .vcenter}


```{r, echo = FALSE, out.width = '85%', fig.retina = NULL}
knitr::include_graphics("./img/Derivative_GIF.gif")
```



## Interpretation

The first derivative tells you about:

- the rate of change
- the slope of the tangent at $x$
- describes the best __linear approximation__ of $f$ near $x$
    - $f'(x) > 0 (< 0)$ then $f$ is increasing (decreasing) around $x$
    - $f'(x) = 0$ then $f(x)$ is a critical point and $x$ is a critical value (the behaviour of $f$ is unclear)

## Second derivative

The second derivative of $f$, denoted $f''$ is the first derivative of $f'$ ($f'' = (f')'$).
The second derivative tells us about:

- describes the change in the rate of change of $f$.
- describes the curvature of $f$ near $x$ (best __quadratic approximation__):
    - $f''(x) > 0 (<0)$ then $f$ is curving upwards (downwards) around $x$
    - $f''(x) = 0$ then the behaviour of $f$ is unclear
    - what if $f'(x^*) = 0$ and $f''(x^*) > 0 (or < 0)?$

## Interpretation of derivatives in physics

If $f(x)$ describes the position of $x$:

- $f'(x)$ describes the __velocity__ of $x$
- $f''(x)$ describes the __acceleration__ of $x$
- Useful for pub quizzes:
    - $f'''(x)$ describes the __jerk__ of $x$
    - $f^{iv}(x)$ describes the __jounce__ of $x$

## Minimum and maximum

When $f'(x) = 0$, $f(x)$ is a critical point:

- if $f''(x) > 0$ then $f(x)$ is a local minimum of $f$
- if $f''(x) < 0$ then $f(x)$ is a local maximum of $f$


## Partial derivative

For multivariate function, such as $f(x, y)$, a partial derivative is the derivative of this function with respect to one variable, with others held constant. It is usually noted $f'_x(x, y)$ (another common notation is $\frac{\partial f}{\partial x}$).

```{r, echo = FALSE, out.width = '50%', fig.retina = NULL, fig.align = 'center'}
knitr::include_graphics("./img/Partial_func_eg.png")
```

## Properties of partial derivative

- One can define the second partial derivative of $f$, noted $\frac{\partial^2}{\partial x_i \partial x_j} f$
- Schwarz' theorem: $\frac{\partial}{\partial x_i \partial x_j} f = \frac{\partial}{\partial x_j \partial x_i} f$ if the second derivatives are continuous



## Gradient

The gradient of a multivariate function $f: \mathbb{R}^n \mapsto \mathbb{R}$  is noted $\nabla f$ and defined as the vector of partial derivatives of $f$:

$$
\nabla f(a) = \begin{bmatrix} \frac{\partial}{\partial x_1} f(a)\\
                     \vdots\\
                     \frac{\partial}{\partial x_n} f(a)
              \end{bmatrix}
$$

$\nabla f$ points in  the direction of the greatest rate of increase of $f$

The magnitude of $\nabla f$ is the slope of the graph in that direction

## Hessian

The Hessian (or Hessian matrix) is a square matrix of the partial second derivatives of a function.

$$
H(a) = \begin{bmatrix}
              \frac{\partial^2}{\partial^2 x_{1}} f(a) & \frac{\partial^2}{\partial x_1\partial x_2} & \cdots & \frac{\partial^2}{\partial x_1\partial x_n}\\
              \frac{\partial^2}{\partial x_2\partial x_1} f(a) & \frac{\partial^2}{\partial^2 x_{2}} f(a) & \cdots & \frac{\partial^2}{\partial x_{2} \partial x_n} f(a)\\
              \vdots & \vdots & \ddots & \vdots\\
              \frac{\partial^2}{\partial x_n \partial x_1} f(a) & \frac{\partial^2}{\partial x_n \partial x_2} & \cdots & \frac{\partial^2}{\partial^2 x_n}
       \end{bmatrix}
$$
If all second derivates are continuous, the Hessian is symmetric (Schwartz' theorem)




## Some example of functions

The constant function takes the form: $f(x) = a$

A simple linear function: $f(x) = ax + b$

## $f(x) = a$

```{r, out.width = '60%', fig.retina = NULL, fig.align = "center"}
library(ggplot2)
data <- data.frame(x = 1:10, fx = rep(5, 10))
ggplot(data, aes(x = x, y = fx)) + geom_line(colour="#FF3465")
```

## $f(x) = ax + b$

```{r, out.width = '60%', fig.retina = NULL, fig.align = "center"}
x <- 1:10
fx <-  5 * x + 2
data <- data.frame(x = x, fx = fx)
ggplot(data, aes(x = x, y = fx)) + geom_line(colour="#FF3465")
```

## $f(x) = a x^2 + b$

```{r, out.width = '60%', fig.retina = NULL, fig.align = "center"}
x <- 1:10
fx <-  2 * x ^ 2 - 10
data <- data.frame(x = x, fx = fx)
ggplot(data, aes(x = x, y = fx)) + geom_line(colour="#FF3465")
```


## Multivariate functions 

```{r, out.width = '60%', fig.retina = NULL, fig.align = "center"}
x <- 1:35; y <- 5:35
x_seq <- seq(min(x), max(x), length.out = 40)
y_seq <- seq(min(y), max(y), length.out = 40)
fx <- outer(x_seq, y_seq, function(x, y){ x ^ 2 - 5 * y + 3})
persp(x = x_seq, y = y_seq, z = fx, theta = -30, phi = 30, 
      xlab = "x", ylab = "y", zlab = "f(x)", col = "blue", expand = 0.8)
```

## RGL

```{r, webgl = TRUE, fig.align = "center"}
#library(rgl)
f <- function(x, y){ x ^ 2 - 5 * y + 3}
persp3d(f)
```

## Derivatives in R

```{r}
D(expression(cos(x)), "x")
D(expression(x^2), "x")
D(D(expression(x^3), "x"), "x")
```

## Derivatives in R (2)

The function `deriv` helps construct a function that will return both $f(x)$ and $f'(x)$:

```{r}
cos_deriv <- deriv(expression(cos(x)), "x", function.arg = TRUE)
cos_deriv(1:5)
```

## Derivatives in R (3)

```{r, fig.align = "center", out.width = "60%"}
cd <- data.frame(x = seq(0, 10, by = 0.01), 
                 cos = as.numeric(cos_deriv(seq(0, 10, by = 0.01))), 
                 cos_prime = as.numeric(attr(cos_deriv(seq(0, 10, by = 0.01)),
                                             "gradient")))
ggplot(data = cd, aes(x = x)) + 
       geom_line(aes(y = cos),color="#FF3489") + 
       geom_line(aes(y = cos_prime),color="#FF7645")
```

## Derivatives in R (4): a better graph

```{r, fig.align = "center", out.width = "60%"}
library(reshape2)
cd2 <- melt(cd, id.vars = "x")
ggplot(data = cd2, aes(x = x, y = value, linetype = variable)) +         geom_line(size = 2,colour="#FF1789")
```

## Hessian

```{r}
library(pracma)
f <- function(x) cos(x[1] + x[2])
x0 <- c(0, 0)
hessian(f, x0)
```

## Hessian (2)


```{r}
f <- function(u) {
    x <- u[1]; y <- u[2]; z <- u[3]
    return(x^3 + y^2 + z^2 +12*x*y + 2*z)
}
x0 <- c(1,1,1)
hessian(f, x0)
```


# {.flexbox .vcenter}

![](img/fin.png)