\documentclass[a4paper]{article}

\usepackage[english]{babel}
\usepackage{amsmath}
\usepackage{amsfonts}
\usepackage{amssymb}
\usepackage{fullpage}
\usepackage{mathtools}
\usepackage{float}
\usepackage{enumerate}
\usepackage{hyperref}
\usepackage{comment}

\setlength{\parindent}{0pt}

\usepackage{amsthm}
\newtheorem{theorem}{Theorem}[section]
\newtheorem*{corollary}{Corollary}
\newtheorem{lemma}{Lemma}[section]
\newtheorem{claim}{Claim}[section]
\newtheorem{axiom}{Axiom}[section]

\theoremstyle{definition}
\newtheorem{definition}{Definition}[section]

\theoremstyle{remark}
\newtheorem{remark}{Remark}[section]
\newtheorem{example}{Example}[section]
\newtheorem{observation}{Observation}[section]

\let\originalleft\left
\let\originalright\right
\renewcommand{\left}{\mathopen{}\mathclose\bgroup\originalleft}
\renewcommand{\right}{\aftergroup\egroup\originalright}
\newcommand{\laplace}[1]{\mathcal{L}\left\{#1\right\}}
\newcommand{\pset}[1]{\mathcal{P}\left(#1\right)}
\newcommand{\pmat}[1]{\begin{pmatrix}#1\end{pmatrix}}
\newcommand{\bmat}[1]{\begin{bmatrix}#1\end{bmatrix}}
\newcommand{\N}{\mathbb{N}}
\newcommand{\R}{\mathbb{R}}
\newcommand{\Z}{\mathbb{Z}}
\newcommand{\Q}{\mathbb{Q}}
\newcommand{\C}{\mathbb{C}}
\renewcommand{\P}{\mathbb{P}}
\newcommand{\inprod}[1]{\left\langle#1\right\rangle}
\DeclareMathOperator{\sgn}{sgn}
\DeclareMathOperator{\val}{val}
\DeclareMathOperator{\ord}{ord}
\DeclareMathOperator{\charac}{char}
\DeclareMathOperator{\im}{im}
\DeclareMathOperator{\id}{id}
\DeclareMathOperator{\intemb}{int}
\DeclareMathOperator{\com}{com}
\DeclareMathOperator{\vecspan}{span}
\DeclareMathOperator{\defequiv}{\mkern4mu :\mkern-8mu \iff}

\newcommand{\ind}[1]{\mathbf{1}_{\left\{#1\right\}}}
\newcommand{\norm}[1]{\left\lVert#1\right\rVert}
\newcommand{\hateq}{\mathrel{\widehat{=}}}

\makeatletter
\newcommand{\axcirc}[1]{\mathbin{\mathpalette\make@circled#1}}
\newcommand{\make@circled}[2]{%
	\ooalign{$\m@th#1\largecirc{#1}$\cr\hidewidth$\m@th#1#2$\hidewidth\cr}%
}
\newcommand{\largecirc}[1]{%
	\vcenter{\hbox{\scalebox{1.5}{$\m@th#1\bigcirc$}}}%
}
\makeatother

\newcommand{\letprob}{Let $(S, \Sigma, P)$ be a probability space. {}}
\let\letprobalways\letprob

% (un)comment to remove useless fluff
\renewcommand{\letprob}{}

\let\tilde\widetilde
\let\hat\widehat

\title{Mathematical Statistics}
\author{Koen Oostveen}

\begin{document}
	\begin{titlepage}
		\maketitle
		\thispagestyle{empty}
		\tableofcontents
	\end{titlepage}
	
	\part{Mathematical Statistics I}
	\section{Probability Theory}
	\subsection{Probability spaces}
	\begin{definition}
		A \textbf{sample space} is a nonempty set representing all possible outcomes of an experiment.
	\end{definition}
	\begin{axiom}
		Let $S$ be a sample space. Let $\Sigma \subseteq \pset{S}$ be a $\sigma$-algebra on $S$. Let $P : \Sigma \to [0, 1]$. A triple $(S, \Sigma, P)$ is then called a \textbf{probability space} (the $\Sigma$ is commonly omitted) if the following axioms hold:
		\begin{enumerate}
			\item $\forall E \in \Sigma : P(E) \geq 0$
			\item $P(S) = 1$
			\item Let $(E_i)_{i \in \N^*}$ be a sequence in $\Sigma$. Suppose that for all $i, j \in \N^*$ where $i \neq j$ we have $E_i \cap E_j = \varnothing$. Then
			\[ P\left( \bigcup_{i \in \N^*} E_i\right) = \sum_{i \in \N^*} P(E_i) \]
		\end{enumerate}
	\end{axiom}
	\begin{definition}
		\letprobalways Let $E \in \Sigma$ ($E$ is called an \textbf{event}). We then define the set $E^C$, the so-called \textbf{complement} of $E$ in $S$, by
		\[ E^C := S \setminus E \]
	\end{definition}
	\begin{definition}
		\letprobalways Let $E, F \in \Sigma$ be events. We then define the \textbf{conditional probability} of $E$ \textbf{given} $F$ $P(E|F)$ by
		\[ P(E|F) := \frac{P(EF)}{P(F)} \]
	\end{definition}
	\begin{definition}
		\letprobalways Let $E, F \in \Sigma$. We say that $E$ and $F$ are \textbf{independent} if
		\[ P(EF) = P(E) P(F) \]
	\end{definition}
	\begin{definition}
		\letprobalways Let $E_1, \ldots, E_n \in \Sigma$, $n \in \N^*$. We say that these events are \textbf{pairwise independent} if for all $i, j = 1, \ldots, n$, where $i \neq j$, that $E_i$ and $E_j$ are independent.
		
		We say that these events are \textbf{independent} if
		\[ P\left( \bigcap_{i=1}^n E_i \right) = \prod_{i=1}^n P(E_i) \]
	\end{definition}
	
	\newpage
	\subsection{Random variables}
	\subsubsection{Definitions}
	\begin{definition}
		\letprobalways A map $X : S \to \R$ is called a \textbf{random variable}. We define the set $S_X := \im(X) \subseteq \R$. Let $\Sigma_X \subseteq \pset{S_X}$. Define the map $P_X : \Sigma_X \to \R$ by
		\[ P_X(D) := P(\left \{ s \in S \mid X(S) \in D \right\}) \]
		Then also the triple $(S_X, \Sigma_X, P_X)$ is a probability space.
	\end{definition}
	\begin{definition}
		\letprobalways Let $X$ be a random variable in said space. We define a map $F_X : \R \to \R$ by
		\[ a \mapsto P(X \leq a) := P_X((-\infty, a]) \]
		$F_X$ is called the \textbf{cumulative distribution function} of $X$.
	\end{definition}
	\begin{definition}
		\letprobalways Let $X$ be a random variable in said space. Suppose $X$ is discrete, that is, $S_X$ is finite. Then we define the \textbf{probability mass function} of $X$ $p_x : \R \to \R$ such that
		\[ p_x(a) := P(X=a) \]
	\end{definition}
	\subsubsection{Discrete random variables}
	\begin{definition}
		\letprob Let $X$ be a discrete random variable. Let $p \in [0, 1]$. Suppose $S_X = \{0, 1\}$ and $p(1) = p$, $p(0) = 1-p$. Then $X$ is said to be \textbf{Bernoulli distributed}.
	\end{definition}
	\begin{definition}
		\letprob Let $X$ be a discrete random variable. Let $p \in [0, 1]$, $n \in \N^*$ `the amount of trials'. Suppose $S_X = \{0, \ldots, n\}$ and
		\[ p(i) = \binom{n}{i} p^i (1-p)^{n-i} \]
		Then $X$ is said to be \textbf{Binomially distributed}, and we write $X \sim B(n, p)$.
	\end{definition}
	\begin{definition}
		\letprob Let $X$ be a discrete random variable. Let $p \in [0, 1]$. Suppose $S_X = \N^*$ ($X$ is the `amount of experiments until 1 success') and
		\[ p(i) = p(1-p)^{i-1} \]
		Then $X$ is said to be \textbf{Geometrically distributed}, and we write $X \sim \operatorname{Geo}(p)$. Some properties are:
		\begin{itemize}
			\item $P(X > i) = (1-p)^i$
			\item $P(X > j + i \mid X > j) = P(X > i) = (1-p)^i$
		\end{itemize}
	\end{definition}
	\begin{definition}
		\letprob Let $X$ be a discrete random variable. Let $\lambda > 0$ (`the amount of expected successful events'). Suppose $S_X = \N$ ($X$ is the `amount of experiments that are successful') and
		\[ p(i) = \exp(-\lambda) \frac{\lambda^i}{i!} \]
		Then $X$ is said to be \textbf{Poisson distributed}. 
		
		Suppose that $X \sim B(n, p)$. Let $\lambda := np$. Then
		\[ \lim_{n \to \infty} P(x = i) = p(i) \]
		where $p$ refers to the probability mass function of a Poisson distributed variable. That is, if you do arbitrarily many experiments, a binomially distributed variable behaves like a Poisson variable.
	\end{definition}
	
	\subsubsection{Continuous random variables}
	\begin{definition}
		\letprob Let $X$ be a \textbf{continuous} random variable, that is, $S_X$ is not finite (and maybe also not countable). Let $f : \R \to \R$ be an (improperly) integrable function. Let $B \subseteq S_X$ be a closed interval. If we have the properties that:
		\begin{itemize}
			\item $P_X(B) = P(X \in B) = \int_B f(x) \, dx$
			\item $\forall x \in \R : f(x) \geq 0$
			\item $\int_\R f(x) \, dx = 1$
		\end{itemize}
		Then $f$ is said to be a \textbf{probability density function} of $X$. If $B = [a, b]$ for some $a < b$ real numbers, it follows that
		\begin{itemize}
			\item $P(a \leq X \leq b) = \int_a^b f(x) \, dx$
			\item $F(b) := P(x \leq b) = \int_{-\infty}^b f(x) \, dx$
			\item $\forall x \in \R : f(x) = F'(x)$ (if $f$ is continuous at $x$, which we can assume most of the time)
		\end{itemize}
	\end{definition}
	\begin{definition}
		\letprob Let $X$ be a continuous random variable and suppose it has probability density function $f$. The following will be a list of potential identities for $f$ and associated `labels'. Let $x \in \R$.
		\begin{itemize}
			\item Let $\alpha < \beta$ be real numbers. We say that $X$ is \textbf{uniformly distributed} ($X \sim U(\alpha, \beta)$) if
			\[ f(x) = \begin{cases*}
				\frac{1}{\beta - \alpha} & if $\alpha < x < \beta$ \\
				0 & otherwise
			\end{cases*} \]
			\item Let $\lambda > 0$. We say that $X$ is \textbf{exponentially distributed} if for $x \geq 0$
			\[ f(x) = \lambda \exp(-\lambda x) \]
			where we make the convention that if the density is unspecified somewhere, it is $0$. Observe that for some $b > 0$ we have
			\[ \int_{-\infty}^b f(x) \, dx = \lambda \int_0^b \exp(-\lambda x) = \exp(-\lambda x)\biggr\rvert_b^0 = 1 - \exp(-\lambda b) \]
			Hence indeed
			\[ \int_\R f(x) \, dx = \lim_{b \to \infty} 1 - \exp(-\lambda b) = 1 \] 
			\item Let $\sigma > 0$. Let $\mu \in \R$. We say that $X$ is \textbf{Normally distributed} ($X \sim N(\mu, \sigma^2)$) if
			\[ f(x) = \frac{1}{\sigma \sqrt{2 \pi}} \exp\left (-\frac{(x - \mu)^2}{2 \sigma^2}\right) \]
			We also define the \textbf{standard normal distribution} $\phi$ according to $N(0, 1)$, e.g.
			\[ \phi(x) := \frac{1}{\sqrt{2\pi}} \exp\left (-\frac{x^2}{2}\right) \]
		\end{itemize}
	\end{definition}
	
	\subsubsection{Expectation}
	\begin{definition}
		\letprob Let $X$ be a discrete random variable. Let $S_X \subset \R$ be the range of $X$, which is assumed to be finite (by assumption of a discrete random variable). We define the \textbf{expected value} of $X$, $EX$, by
		\[ EX := \sum_{x \in S_X} x \underbrace{P(X = x)}_{p(x)} \]
	\end{definition}
	\begin{definition}
		\letprob Let $X$ be a continuous random variable. Let $f : \R \to \R$ be the probability density function of $X$, which is assumed to exist. We define the \textbf{expected value} of $X$, $EX$, by
		\[ EX := \int_\R x f(x) \, dx \]
	\end{definition}
	
	\subsubsection{Functions of random variables}
	\begin{definition}
		\letprob Let $X$ be a random variable. Let $S_X \subseteq \R$ be its range. Let $S_Y \subseteq \R$ be a set. Let $g : S_X \to S_Y$ be a surjective function. Let $Y : S \to \R$ be a random variable, such that $Y(s) := g(X(s))$. We say that such random variable $Y$ is the \textbf{mapped random variable} of $X$ under $g$. We write $Y = g(X)$ (by slight abuse of notation).
	\end{definition}
	\begin{definition}
		\letprob Let $X$ be a random variable. We say that for some $n \in \N$ that its $n$-th \textbf{moment} is the expected value of $X^n$, that is, the expected value of a random variable $S \to \R$ such that $s \mapsto X(s)^n$.
		
		We say that for some $n \in \N$ that its $n$-th \textbf{central moment} is the value
		\[ E((X - E(X))^n) \]
		
		We say that the \textbf{variance} of $X$, denoted by $\operatorname{var}(X)$, is the value
		\[ \operatorname{var}(X) := E((X - E(X))^2) \]
		That is, the variance of $X$ is its $2$-nd central moment.
		
		We say that the \textbf{standard deviation} $\sigma_X$ of $X$ is
		\[ \sigma_X := \sqrt{\operatorname{var}(X)} \]
	\end{definition}
	\begin{theorem}
		\letprob Let $X$ be a random variable. Then
		\[ \operatorname{var}(X) = E(X^2) - (E(X))^2 \]
	\end{theorem}
	
	\newpage
	
	\subsection{Joint distributions}
	\subsubsection{Joint distribution functions}
	\begin{definition}
		\letprobalways Let $X$ and $Y$ be discrete random variables in said space. We introduce the following notation: for $i \in S_X$ and $j \in S_Y$,
		\[ P(X = i, Y = j) := P(\{ s \in S \mid X(s) = i \land Y(s) = j \}) \]
		(likewise for multiple variables).
		The function $p : S_X \times S_Y \to [0, 1]$ such that $i, j \mapsto P(X = i, Y = j)$ is called the \textbf{joint probability mass function} of $X$ and $Y$.
	\end{definition}
	\begin{theorem}
		\letprobalways Let $X$ and $Y$ be discrete random variables in said space. Let $p : S_X \times S_Y$ be the joint probability mass function of $X$ and $Y$. Then for all $i \in S_X$
		\[ P(X = i) = \sum_{j \in S_Y} P(X = i, Y = j) \]
		and for all $j \in S_Y$
		\[ P(Y = j) = \sum_{i \in S_X} P(X = i, Y = j) \]
	\end{theorem}
	\begin{definition}
		\letprobalways Let $X$ and $Y$ be continuous random variables in said space. A function $f : \R^2 \to \R$ is said to be the \textbf{joint probability density function} of $X$ and $Y$ if for all (measurable) sets $A \subseteq S_X$ and $B \subseteq S_Y$, we have
		\[ P(X \in A, Y \in B) = \int_B\int_A f(x, y) \, dx \, dy \]
		and for all $x, y \in \R^2$, $f(x, y) \geq 0$.
	\end{definition}
	
	\subsubsection{Independent random variables}
	\begin{definition}
		\letprob Let $X$ and $Y$ be discrete random variables. We say that $X$ and $Y$ are \textbf{independent} if for all $x \in S_X$ and $y \in S_Y$ we have
		\[ P(X = x, Y = y) = P(X = x) P(Y = y) \]
	\end{definition}
	\begin{definition}
		\letprob Let $X$ and $Y$ be random variables. We say that $X$ and $Y$ are \textbf{(mutually) independent} if for all $a, b \in \R$ we have
		\[ P(X \leq a, Y \leq b) = P(X \leq a)P(Y \leq b) \]
	\end{definition}
	\begin{remark}
		If $X, Y$ are discrete random variables then the above definition is equivalent to: for all $x \in S_X$, $y \in S_Y$, we have
		\[ P(X = x, Y = y) = P(X = x) P(Y = y) \]
		and if $X, Y$ are continuous, then: if $f : \R^2 \to \R$ is the density function of $X$ and $Y$, then for all $x, y$ we have
		\[ f(x, y) = f_X(x) f_Y(y) \]
		where $f_X$ and $f_Y$ are the density functions of $X$ and $Y$, respectively.
	\end{remark}
	\begin{theorem}
		\letprob Let $X$ and $Y$ be random variables. If $X$ and $Y$ are independent, then for functions $g : S_X \to \R$ and $h : S_Y \to \R$ we have
		\[ E(g(X)h(Y)) = E(g(X))E(h(Y)) \]
	\end{theorem}
	\begin{corollary}
		If $X, Y$ are independent then
		\[ E(XY) = E(X)E(Y) \]
		The converse is not necessarily true.
	\end{corollary}
	
	\subsection{Limit theorems}
	\begin{theorem}
		\letprobalways Let $X$ be a random variable in said space. If $X \geq 0$ for all outcomes in $S$, then for all $a > 0$ we have
		\[ P(X \geq a) \leq \frac{1}{a} EX \]
	\end{theorem}
	\begin{theorem}
		\letprob Let $X$ be a random variable. If $EX = \mu$ and $\operatorname{var}(X) = \sigma^2$, then for all $k > 0$ we have
		\[ P(|X - \mu| \geq k) \leq \frac{\sigma^2}{k^2} \]
	\end{theorem}
	\begin{proof}
		We already know from the previous theorem that
		\[ P(|X - \mu| \geq k) = P((X - \mu)^2 \geq k^2) \leq \frac{E(X - \mu)^2}{k^2} = \frac{\sigma^2}{k^2} \]
	\end{proof}
	\begin{theorem}
		\letprobalways Let $(X_k)_{k \in \N^*}$ be a sequence of random variables in that space, independent and identically distributed (iid). Suppose for all $k \in \N^*$ we have $E X_k = \mu$ for some $\mu \in \R$. We define a new sequence of random variables, for $n \in \N^*$,
		\[\overline{X}_n := \frac{1}{n} \sum_{k=1}^n X_k\] Then, with probability (wp) 1
		\[ \lim_{n \to \infty} \overline{X}_n = \mu \]
		(whatever that means, to be honest).\\
		Instead of proving this result, we will consider a more tangible / defined result: for all $\varepsilon > 0$ we have that
		\[ \lim_{n \to \infty} P(|\overline{X}_n - \mu| < \varepsilon) = 1 \]
	\end{theorem}
	\begin{theorem}[Central limit]
		\letprob Let $(X_k)_{k \in \N^*}$ be a sequence of random variables, independent and identically distributed (iid). Suppose for all $k \in \N^*$ we have $E X_k = \mu$ for some $\mu \in \R$ and $\operatorname{var}(X_k) = \sigma^2$ for $\sigma > 0$. Then
		\[ P\left( \frac{\overline{X}_n - \mu}{\sigma/\sqrt{n}} \leq a \right) \xrightarrow{n \to \infty} \Phi(a) \]
	\end{theorem}
	
	\begin{definition}
		\letprob Let $(X_k)_{k \in \N^*}$ be a sequence of random variables. We say that $(X_k)$ converges to a random variable $X$ (distributionally) if
		\[ \forall \varepsilon > 0 : \lim_{k \to \infty} P(|X_k - X| \geq \varepsilon) = 0 \]
		or equivalently
		\[ \forall \varepsilon > 0 : \lim_{k \to \infty} P(|X_k - X| < \varepsilon) = 1 \]
	\end{definition}
	
	\begin{theorem}
		Distributional limits are unique (with probability 1).
	\end{theorem}
	\begin{proof}
		Let $(X_k)$ be a sequence of random variables and let $X$ and $\tilde X$ be distributional limits of $(X_k)$. It suffices to show that $P(X = Y) = 1$, that is, $P(\{s \in S \mid X(s) = Y(s)\}) = 1$. We have
		\[ |X - \tilde X| = |X - X_k + X_k - \tilde X| \leq |X - X_k| + |X_k - \tilde X| \]
		Let $x > 0$, therefore
		\[ P(|X - \tilde X| \geq 2x) \leq P(|X_n - X| \geq x) + P(|X_n - \tilde X| \geq x) \]
		By the squeeze theorem and taking limits we find $P(|X - \tilde X| \geq 2x) = 0$, which implies the goal.
	\end{proof}
	
	\newpage
	
	\subsection{$\chi^2$ distribution}
	\begin{definition}
		A $\chi^2$ distribution with $n \in \N^*$ \textbf{degrees of freedom} is the random variable $X$ where $X = X_1^2 + \ldots + X_n^2$, where $X_1, \ldots, X_n$ are i.i.d. and $X_i \sim N(0, 1)$ for all $i$. We write $X \sim \chi^2(n)$.
	\end{definition}
	\begin{remark}
		Observe that for such $X \sim \chi^2(n)$ we have
		\[ EX = E\left[ X_1^2 + \ldots + X_n^2 \right] = EX_1^2 + \ldots + EX_n^2 = n \]
		Similarly, for any $i$:
		\begin{align*}
			\operatorname{Var}(X_i^2) &= EX_i^4 - (EX_i^2)^2 \\
			&= \frac{1}{\sqrt{2 \pi}} \int_\R x^4 \exp\left(-\frac{1}{2} x^2\right) \, dx - 1 \\
			&= \frac{1}{\sqrt{2 \pi}} \left( \int_\R 3x^2 \exp\left(-\frac{1}{2} x^2\right) \, dx \right) - 1 \\
			&= 3 \left( \frac{1}{\sqrt{2 \pi}} \int_\R x^2 \exp\left(-\frac{1}{2} x^2\right) \, dx \right) - 1 \\
			&= 3 - 1 = 2
		\end{align*}
		Therefore
		\[ EX = \operatorname{Var}(X_1^2 + \ldots + X_n^2) = 2n \]
		It turns out this distribution is already known, namely $X \sim \operatorname{Gamma}(n/2, 2)$. Therefore, for even $n = 2(k + 1)$, $k \in \N$, we get the following `nice' formula for the density
		\[ f_X(x) := \begin{cases*}
			\dfrac{\exp\left( -\frac{x}{2} \right) \left(\frac{x}{2}\right)^k}{k!} & if $x > 0$ \\
			0 & otherwise
		\end{cases*} \]
	\end{remark}
	
	\newpage
	
	\section{Descriptive statistics}
	\subsection{Order statistics}
	\begin{definition}
		\letprobalways Let $X_1, \ldots, X_n$ be random variables in $S$. We define now new random variables, $X_{(1)}, \ldots, X_{(n)}$ in $S$, with the following property: using the natural ordering of the real numbers, we can, given some outcome $s \in S$, find numbers $X_{(i)}(s)$ for $i = 1, \ldots, n$ such that
		\[ X_{(1)}(s) \leq \ldots \leq X_{(n)}(s) \]
		Notice that there is guaranteed to be a permutation, that is, bijective $p : \{1, \ldots, n\} \to \{1, \ldots, n\}$, such that
		\[ X_{p(1)}(s) \leq \ldots \leq X_{p(n)}(s) \]
		however, such $p$ might not be unique if there are duplicate values. However, obviously, independent on valid choice of $p$, this defines a unique sequence $X_{p(1)}(s), \ldots, X_{p(n)}(s)$.
		We then define the random variables each $X_{(i)}$ pointwise by choosing a valid $p$ for every $s \in S$ and setting $X_{(i)}(s) := X_{p(i)}$.\\
		The now defined random variables are called the \textbf{order statistics} of $X_1, \ldots, X_n$.
	\end{definition}
	\begin{remark}
		Observe that these order statistics are generally very dependent. One can see this by choosing two independent $X_1, X_2 \sim B(p)$ for some $0 < p < 1$. We can compute the distributions of the order statistics. There are only the following possibilities for $(X_1, X_2)$: $(0, 0)$, $(0, 1)$, $(1, 0)$, $(1, 1)$. Therefore:
		\[ P(X_{(1)} \leq x) = P(\min\{X_1, X_2\} \leq x) = 1 - P(\min\{X_1, X_2\} > x) = 1 - P(X_1 > x)P(X_2 > x) \]
		\[ P(X_{(2)} \leq x) = P(\max\{X_1, X_2\} \leq x) = P(X_1 \leq x) P(X_2 \leq x) \]
		We can say something about independence:
		\[ P(X_{(1)} \leq 0, X_{(2)} \leq 0) = (1-p)^2 \]
		while
		\[ P(X_{(1)} \leq 0) P(X_{(2)} \leq 0) = (1-p^2)(1-p)^2 \]
		therefore, we get independence only if $1 - p^2 = 1$, which is false.
	\end{remark}
	\begin{definition}
		Let $X_1, \ldots, X_n$ be random variables with order statistics $X_{(1)}, \ldots, X_{(n)}$. The \textbf{median} is then a new random variable $M$ such that
		\begin{itemize}
			\item If $n$ is odd, $M = X_{\left(\frac{n+1}{2}\right)}$
			\item If $n=2k$ (where $k \in \N$) is even, $\displaystyle M = \frac{X_{(k)} + X_{(k + 1)}}{2}$
		\end{itemize}
	\end{definition}
	\subsection{Quantile function}
	\begin{remark}
		Below follows a stupid ill-defined thing, and I'll say why it is ill-defined soon enough.
	\end{remark}
	\begin{definition}
		Let $F : \R \to [0, 1]$ be a continuous c.d.f. of a random variable $X$, and suppose it is bijective. The function $Q := F^{-1}$ is called the \textbf{quantile function} of $X$.
	\end{definition}
	\begin{remark}
		Issue: a continuous bijection $\R \to [0, 1]$ could not exist, if it did then the inverse is also continuous (because we are dealing with $\R$), and then $Q^{-1}([0, 1]) \neq \R$ because $[0, 1]$ is compact and $\R$ is not. Therefore, no such $F$ could exist. Therefore, we adopt the following definition.
	\end{remark}
	\begin{definition}
		Let $F : \R \to [0, 1]$ be a continuous c.d.f. of a random variable $X$. We define a new function, the \textbf{quantile function} $Q : [0, 1] \to \R$ of $X$, by $p \mapsto Q(p)$, where
		\[ Q(p) := \inf \{ x \in \R \mid p \leq \underbrace{F(x)}_{P(X \leq x)} \} \]
	\end{definition}
	\begin{remark}
		This function $Q$ satisfies now being an `inverse' of $F$ in the sense that
		\[ \forall p \in [0, 1] : F(Q(p)) = P(X \leq Q(p)) = p \]
		i.e. $F \circ Q = \id_{[0, 1]}$, which is better known as $Q$ being a \textbf{right inverse} of $F$. Notice that we could have chosen other functions $Q$ that satisfy being a right inverse.
	\end{remark}
	\begin{definition}
		Let $X$ be a random variable with quantile function $Q$. Let $k = 1, \ldots, 100$. The $k$-th \textbf{percentile} is defined to be the number $Q(k/100)$, i.e. the number $q \in \R$ such that $k\%$ of the time, $X \leq q$.
	\end{definition}
	\begin{definition}
		Let $X = (X_1, \ldots, X_n) \in \R^n$ be a \textbf{sample}, i.e. values that corresponding i.i.d. random variables took for a certain outcome $s \in S$. We define the following quantities:
		\begin{itemize}
			\item The \textbf{total} $\displaystyle \sum X := \sum_{i=1}^n X_i$
			\item The \textbf{sample mean} $\displaystyle \bar X := \frac{\sum X}{|X|}$
			\item The \textbf{sample median}, by the value $M(s)$, where $M$ is the random variable corresponding with the median of the original random variables (i.e. it is defined similarly).
			\item The \textbf{sample variance} $s_X^2 := \frac{1}{n-1} \sum_{i=1}^n (X_i - \bar X)^2$.
			\item The \textbf{sample standard deviation} $s_X := \sqrt{s_X^2}$.
			\item The \textbf{interquartile range} $Q_3 - Q_1 = Q(3/4) - Q(1/4)$.
		\end{itemize}
	\end{definition}
	\begin{remark}
		Why is the sample variance defined that way? Well, if we let
		\[ S^2 := \frac{1}{n} \sum_{i=1}^n (X_i - \bar X)^2 \]
		we can compute that
		\[ ES^2 = \frac{n-1}{n} \sigma^2 \]
		Therefore, the better approximation that we can make of $\sigma^2$ is
		\[ s^2 := \frac{n}{n-1} S^2 = \frac{1}{n-1} \sum_{i=1}^n (X_i - \bar X)^2 \]
		to remove the `bias'. By the way, this definition of $s^2$ is also a random variable, but if you close your eyes it becomes a realization when you confuse the i.i.d. $X_1, \ldots, X_n$ with a sample of it.
	\end{remark}
	\begin{definition}
		Let $X$ be a random variable. Let $k \in \N^*$. We say that the \textbf{k-th corrected central moment} of $X$ is the quantity
		\[ \frac{E(X - EX)^k}{\sigma^k} \]
		where $\sigma := \sqrt{E(X - EX)^2}$. We define the additional quantities:
		\begin{itemize}
			\item The \textbf{skewness} of $X$ is this moment when $k=3$
			\item the \textbf{kurtosis} $K_X$ of $X$ is this moment when $k=4$.
			\item The \textbf{excess kurtosis} of $X$ is $K_X - 3$.
		\end{itemize}
	\end{definition}
	\begin{remark}
		Why such a weird last definition? Well, a normally distributed $X \sim N(\mu, \sigma^2)$ has a kurtosis of $3$, kurtosis in a way measures how much `tail' a distribution has, so if it has a higher kurtosis than the normal distribution, you can expect a larger number of statistical outliers than with the normal distribution.
	\end{remark}
	\section{Estimators}
	\begin{definition}
		Let $X_1, \ldots, X_n$ be a random sample (i.i.d.), then a (measurable) $T : \R^n \to \R$ is called a \textbf{statistic}, which is also often conflated with the mapped random variable $T(X_1, \ldots, X_n)$, also called a \textbf{statistic}.
	\end{definition}
	\begin{definition}
		Let $\Theta$ be any set, referred to as the parameter space. Actually, in this course we only care about $\Theta = \R$. An \textbf{estimator} $T$ is then a statistic on any finite random sample.
	\end{definition}
	\begin{definition}
		Let $\hat \theta$ be an estimator of $\theta$. The \textbf{bias} of $\hat \theta$, $\operatorname{Bias}(\hat \theta)$, is defined by
		\[ \operatorname{Bias}(\hat \theta) := E \hat \theta - \theta \]
	\end{definition}
	\begin{definition}
		Let $\theta \in \R$ be a parameter to be estimated. An estimator $T$ is said to be \textbf{unbiased} if $ET = \theta$ (in the above sense where $T$ is taken to be a random variable).
	\end{definition}
	\begin{example}
		$s^2$ is an unbiased estimator of $\sigma$.
	\end{example}
	\begin{proof}
		We established that $Es^2 = \frac{n}{n-1} ES^2 = \sigma^2$ earlier.
	\end{proof}
	\begin{remark}
		Even though $s^2$ is an unbiased estimator of $\sigma^2$, if $\mu$ is known the following is also an unbiased estimator of $\sigma^2$:
		\[ \hat \sigma^2 := \frac{1}{n} \sum_{i=1}^n (X_i - \mu)^2 \]
	\end{remark}
	\begin{example}
		The estimator $\hat \mu := \frac{1}{n} \sum_{i=1}^n X_i$ is an unbiased estimator for $EX$.
	\end{example}
	\begin{proof}
		$\displaystyle E\left[ \frac{1}{n} \sum_{i=1}^n X_i \right] = \frac{1}{n} nEX = EX$
	\end{proof}
	\begin{example}
		We also have that $\operatorname{Var}(\hat \mu) = \sigma^2/n$, i.e. if $n$ increases, the inaccuracy of the random variable $\hat \mu$ becomes less and less in a sense. We can put more rigor to this by observing that
		\[ P(|\hat \mu - \mu| \geq c) = P(|\hat \mu - E\hat \mu| \geq c) \leq \frac{\operatorname{Var}(\hat \mu)}{c^2} = \frac{\sigma^2}{nc^2} \]
		Therefore, in particular for normally distributed variables, we can ask for the probability of the sample mean being more than a stdev away from the real mean
		\[ P(|\hat \mu - \mu| \geq \sigma) \leq \frac{1}{n} \]
		This also means that indeed $\hat \mu_n \to \mu$ as $n \to \infty$ in probability.
	\end{example}
	\begin{example}
		Let $X \sim B(n, p)$, and suppose $n$ is known. An unbiased estimator of $p$ is $\hat p(X) := \frac{X}{n}$.
	\end{example}
	\begin{proof}
		$\displaystyle E\hat p = \frac{EX}{n} = \frac{np}{n} = p$.
	\end{proof}
	\begin{remark}
		Unbiased estimators are not necessarily \textbf{good} estimators. Hence the following definition in order to measure how `good' an estimator is (analogous to the second central moment).
	\end{remark}
	\begin{definition}
		Let $\theta \in \Theta \subseteq \R$ be a parameter to be estimated and let $\hat \theta : \R^n \to \Theta$ be an estimator. The \textbf{mean square error} of $\hat \theta$, $\operatorname{MSE}(\hat \theta)$, is defined by
		\[  \operatorname{MSE}(\hat \theta) := E[(\hat \theta - \theta)^2] \]
	\end{definition}
	\begin{lemma}
		For all estimators $\theta$, $\displaystyle \operatorname{MSE}(\hat \theta) = \operatorname{Bias}(\hat \theta)^2 + \operatorname{Var}(\hat \theta)$.
	\end{lemma}
	\begin{proof}
		\begin{align*}
			E[(\hat \theta - \theta)^2] &= E[\hat \theta^2 - 2 \hat \theta \cdot \theta + \theta^2] \\
			&= E\hat \theta^2 - 2 \theta (E \hat \theta) + \theta^2 \\
			&= (E\hat \theta)^2 - 2\theta (E\hat \theta) + \theta^2 + E\hat \theta^2 - (E \hat \theta)^2 \\ 
			&= (E\hat \theta - \theta)^2 + E\hat \theta^2 - (E \hat \theta)^2 \\
			&= \operatorname{Bias}(\hat \theta)^2 + \operatorname{Var}(\hat \theta) \qedhere
		\end{align*}
	\end{proof}
	\begin{example}
		Recall the estimator
		\[ \hat \sigma^2 := \frac{1}{n-1} \sum_{i=1}^n (X_i - \bar X)^2 \]
		is unbiased, but we would like to compare its MSE to what you would get for the estimator
		\[ \tilde \sigma^2 := \frac{1}{n} \sum_{i=1}^n (X_i - \mu)^2 \]
		if $\mu$ is known, which is also unbiased.
		\[ \operatorname{MSE}(\hat \sigma^2) = \operatorname{Var}(\hat \sigma^2) = \frac{2\sigma^4}{n-1} \]
		Why? Because $\hat \sigma^2 \frac{n-1}{\sigma^2} \sim \chi^2_{n-1}$ which has expectation $2(n-1)$. Another computation shows
		\[ \operatorname{MSE}(\tilde \sigma^2) = \frac{2\sigma^4}{n} \]
		which is indeed slightly better.
	\end{example}
	
	\newpage
	
	\subsection{Sufficiency}
	
	\begin{definition}
		Let $X_1, \ldots, X_n$ be jointly distributed with some parameter $\theta \in \Theta \subseteq \R$. Let $T~:~\R^n~\to~\R$ be a statistic. $T$ is said to be a \textbf{sufficient statistic} for $\theta$ if the distribution of $X_1, \ldots X_n | \tilde T$ is independent of $\theta$, where $\tilde T$ is the mapped random variable $\tilde T = T(X_1, \ldots, X_n)$.
	\end{definition}
	\begin{remark}
		The above definition is really unclear, since we don't know what independence actually means in this context. Here are two other equivalent formulations:
		\begin{itemize}
			\item Let $f_{X_1, \ldots, X_n} : \Theta \times \R^n \to \R$ be the pdf of $X_1, \ldots, X_n$ with parameter $\theta$. Let $f_{X_1, \ldots X_n | T} : \Theta \times \R^{n+1} \to \R$ be the conditional pdf of $X_1, \ldots, X_n$ given $\tilde T$ corresponding to $\theta$. Then there must exist a function $f : \R^{n+1} \to \R$ such that for all $\theta$ and points $x_1, \ldots, x_n \in \R$ and $t \in \R$ we have
			\[ f_{X_1, \ldots X_n | T}(\theta, x_1, \ldots, x_n, t) = f(x_1, \ldots, x_n, t) \]
			\item $I(\theta, \tilde T) = I(\theta, X_1, \ldots, X_n)$, where $I$ is the Shannon mutual information, which is to say, you learn as much about $\theta$ by knowing $T$, as by knowing all random variables.
		\end{itemize}
	\end{remark}
	\begin{theorem}
		Another equivalent notion is that there exist nonnegative functions $g : \Theta \times \R \to \R$ and $h : \R^n \to \R$ such that for all $\theta \in \Theta$ and points $x_1, \ldots, x_n \in \R$ we have
		\[ f_{X_1, \ldots, X_n}(\theta, x_1, \ldots, x_n) = h(x_1, \ldots, x_n) g(\theta, T(x_1, \ldots, x_n)) \]
	\end{theorem}
	\begin{proof}
		If you have sufficiency, we may take $g := f_T$,
		where $f_T : \Theta \times \R \to \R$ is the pdf of $T$ depending on $\theta$, and we take
		\[ h(x_1, \ldots, x_n) := f(x_1, \ldots, x_n, T(x_1, \ldots, x_n)) \]
		where $f$ is such that
		\[ f_{X_1, \ldots X_n | T}(\theta, x_1, \ldots, x_n, t) = f(x_1, \ldots, x_n, t) \]
		which clearly recovers the desired result. The other direction is WIP.
	\end{proof}
	\begin{lemma}
		Let $\hat \theta$ be an estimator of $\theta \in \Theta$. Let $T$ be a sufficient statistic for $\theta$. We define the Rao-Blackwell estimator, $\tilde \theta$ of $\hat \theta$, by $\tilde \theta := E[\hat \theta|T]$. Then
		\begin{itemize}
			\item $\displaystyle E[\hat \theta] = E[\tilde \theta]$
			\item $\displaystyle \operatorname{Var}(\tilde \theta) \leq \operatorname{Var}(\hat \theta)$
		\end{itemize}
		i.e. the Rao-Blackwell estimator is never worse than the original estimator, as well as unbiased when the original one is too.
	\end{lemma}
	\begin{definition}
		An estimator $\hat \theta$ of $\theta$ is called \textbf{minimum variance unbiased (MVU)} if it is unbiased and minimizes the variance, in the sense that for any other estimator $\tilde \theta$, $\operatorname{Var}(\hat \theta) \leq \operatorname{Var}(\tilde \theta)$.
	\end{definition}
	\newpage
	\subsection{Common methods to finding estimators}
	There are two main methods, we will outline both briefly.
	\subsubsection{Method of moments}
	Typically, if you are trying to estimate a parameter $\theta$, you can recover $\theta$ through the moments of your distribution. For example, if we have $X \sim N(\mu, \sigma^2)$, we have that $EX = \mu$, so a good first attempt at estimating $\mu$ is $\hat \theta := \bar X$. More generally, if there exists a function $g : \R^n \to \R$ such that $\theta = g(EX, \ldots, EX^n)$, we define the sample moments
	\[ M_i := \frac{1}{n} \sum_{j=1}^n X_j^i \]
	and the estimator $\hat \theta := g(M_1, \ldots, M_n)$.
	\subsection{Maximum likelihood estimation (MLE)}
	Let $L(\cdot|X_1, \ldots, X_n) : \Theta \to [0, \infty)$ be either the likelihood or log-likelihood function. Clearly, a good estimator $\hat \theta$ satisfies that $L(\hat \theta|X_1, \ldots, X_n)$ is large, because that would indicate that somehow the likelihood of $\hat \theta$ being the parameter that was used to sample $X_1, \ldots, X_n$ being large. As an example consider again an i.i.d. sample $X_1, \ldots, X_n \sim N(\mu, 1)$. We compute the density first:
	\begin{align*}
		f_{X_1, \ldots, X_n}(x_1, \ldots, x_n) &= \prod_{i=1}^n f(x_i) \\
		&= \prod_{i=1}^n \frac{1}{\sqrt{2\pi}} \exp\left(- \frac{(x_i - \mu)^2}{2} \right) \\
		&= (2\pi)^{-n/2} \exp\left(- \sum_{i=1}^n \frac{(x_i - \mu)^2}{2} \right)
	\end{align*}
	The log-likelihood is
	\begin{align*}
		L(\theta|x_1, \ldots, x_n) &= \log \left( (2\pi)^{-n/2} \exp\left(- \sum_{i=1}^n \frac{(x_i - \mu)^2}{2} \right)\right) \\
		&= \frac{-n}{2} \log(2\pi) - \frac{1}{2} \sum_{i=1}^n (\mu - x_i)^2
	\end{align*}
	We differentiate
	\begin{align*}
		\frac{\partial}{\partial \mu} L(\theta|x_1, \ldots, x_n) &= - \frac{1}{2} \sum_{i=1}^n \frac{\partial}{\partial \mu} (\mu - x_i)^2 \\
		&= \sum_{i=1}^n (x_i - \mu) = \sum_{i=1}^n x_i - n\mu
	\end{align*}
	Therefore $L$ is maximized in $\mu$ if $\mu = \bar X$, which is the same estimator we already got before.
	\newpage
	\subsection{Completeness}
	\begin{definition}
		A statistic $T$ is called \textbf{complete} if for all measurable functions $h$ and parameters $\theta$, we have that $E_\theta(h(T)) = 0$ implies $h = 0$ a.e. In this course we assume that all functions are measurable and that $h = 0$ as functions.
	\end{definition}
	\begin{example}
		Let $X_1, \ldots, X_n \sim B(p)$, $p \in (0, 1)$. Let $T := \frac{1}{n} \sum_{i=1}^n X_i$. Then $ET = p$, so it is an unbiased estimator for $p$. But is it sufficient and complete? We compute
		\begin{align*}
			P\left(X_1 = x_1, \ldots, X_n = x_n\right) &= \prod_{i=1}^n P(X_i = x_i) \\
			&= \prod_{i=1}^n \begin{cases*}
				p & if $x_i=1$ \\
				1 - p & if $x_i=0$
			\end{cases*} \\
			&= \prod_{i=1}^n p^{x_i} (1-p)^{1-x_i} \\
			&= p^{\sum_{i=1}^nx_i} (1-p)^{n-\sum_{i=1}^nx_i}
		\end{align*}
		so sufficiency is satisfied. What about completeness? Observe that $nT \sim B(n, p)$. Let $h$ be a function. We compute
		\begin{align*}
			E(h(T)) &= \sum_{t \in S_T} h(t) P(T=t) \\
			&= \sum_{t \in S_T} h(t) P(nT=nt) \\
			&= \sum_{t \in S_T} h(t) \binom{n}{nt} p^{nt} (1-p)^{n - nt} \\
			&= (1-p)^{n} \sum_{t \in S_T} h(t) \binom{n}{nt} \left(\frac{p}{1-p}\right)^{nt} = 0
		\end{align*}
		implies clearly that $h = 0$, since this should hold for all $p \in (0, 1)$.
	\end{example}
	\begin{lemma}
		The Rao-Blackwell estimator of a complete statistic is unique (independent of choice of original estimator), provided unbiased estimators of $\theta$ exist.
	\end{lemma}
	\begin{proof}
		Let $T$ be a complete statistic, let $\hat \theta, \hat \theta'$ be unbiased estimators. Let $g : S_T \to \R$ be defined by $t \mapsto E[\hat \theta|T=t] - E[\hat \theta'|T=t]$. Clearly it suffices to show that $g = 0$, in which case $0 = g(T) = \tilde \theta - \tilde \theta'$ which implies $\tilde \theta = \tilde \theta'$, which shows uniqueness. Due to completeness, it suffices to show that $E(g(T)) = 0$, which we have due to unbiasedness of $\hat \theta, \hat \theta'$.
	\end{proof}
	\begin{theorem}
		Let $T$ be a complete and sufficient statistic for $\theta \in \Theta$ on a random sample $X_1, \ldots, X_n$ i.i.d. that has an unbiased estimator. Then the Rao-Blackwell estimator is the unique MVU estimator of $\theta$.
	\end{theorem}
	\begin{proof}
		Due to the lemma it suffices to show that the Rao-Blackwell estimator $\tilde \theta$ is MVU. Let $\hat \theta$ be any unbiased estimator. Then clearly
		\[ \operatorname{Var}(\tilde \theta) = \operatorname{Var}(E[\hat \theta|T]) \leq \operatorname{Var}(\hat \theta) \]
		which finishes the proof.
	\end{proof}
	\newpage
	\subsection{Cramér-Rao lower bound}
	\begin{definition}
		Stupid definitions incoming,
		\begin{itemize}
			\item The \textbf{likelihood function} maps $\theta \in \Theta$ to the random variable $f_X(X|\theta)$, where $x \mapsto f_X(x|\theta)$ is the density of a random variable $X$.
			\item The \textbf{log-likelihood function} maps $\theta \in \Theta$ to $\log f_X(X|\theta) =: L(\theta|X)$.
			\item The \textbf{score function} is the derivative of the log-likelihood function.
		\end{itemize}
	\end{definition}
	\begin{remark}
		How is differentiation defined on functions that map real numbers to random variables? It turns out that we would need a metric on random variables, which arguably we could make, for example we can define $|X - Y| := \sup_{s \in S} |X(s) - Y(s)|$. But in this specific scenario it is not necessary. We simply fix an outcome $s \in S$ and thus an $x = X(s) \in S_X$, then differentiate the function $\theta \mapsto f(x, \theta)$, we then basically get a function $x \mapsto (\theta \mapsto L'(\theta|x))$ which we turn into a random variable $\theta \mapsto L'(\theta|X)$ as such.
	\end{remark}
	\begin{lemma}
		The score function has mean $0$, given dominated convergence / interchangeability of integration and (partial) differentiation.
	\end{lemma}
	\begin{proof}
		\begin{align*}
			E(\ldots) &= E\left( \frac{f'(X|\theta)}{f(X|\theta)} \right) \\
			&= \int_\R  \frac{f'(x|\theta)}{f(x|\theta)} f(x|\theta) \, dx \\
			&= \int_\R f'(x|\theta) \, dx = 0 \qedhere
		\end{align*}
	\end{proof}
	\begin{definition}
		The \textbf{Fischer information} $I_n(\theta)$ of $\theta$ given a random sample $X_1, \ldots, X_n$ is the variance of the score function.
	\end{definition}
	\begin{lemma}
		The Fischer information satisfies (where $s(\theta|X)$ is the score function)
		\[ I_n(\theta) = E(s(\theta|X)^2) \]
	\end{lemma}
	\begin{proof}
		Trivial due to $E(s(\theta|X)) = 0$.
	\end{proof}
	\begin{lemma}
		If the second derivative of $\log L$ exists, and
		\[ \int_\R \frac{\partial^2}{\partial \theta^2} f(x|\theta) \, dx = \frac{\partial^2}{\partial \theta^2} \int_\R f(x|\theta) \, dx \]
		we have that
		\[ I_n(\theta) = -E\left( \frac{\partial^2}{\partial \theta^2} \log f(x|\theta) \right) \]
	\end{lemma}
	\begin{lemma}
		$\displaystyle I_n(\theta) = nI_1(\theta) $
	\end{lemma}
	\begin{proof}
		We proceed by induction. This is very obviously true if $n = 1$. Suppose we have for some $n \in \N$ that $I_n(\theta) = nI_1(\theta)$. We compute
		\begin{align*}
			I_{n+1}(\theta) &= E \left( \left( \frac{\partial}{\partial \theta} \log f(X_1, \ldots, X_{n+1} | \theta) \right)^2 \right) \\
			&= E \left( \left( \frac{\frac{\partial}{\partial \theta} f(X_1, \ldots, X_{n+1} | \theta)}{f(X_1, \ldots, X_{n+1} | \theta)} \right)^2 \right) \\
			&= E \left( \left( \frac{\frac{\partial}{\partial \theta} f(X_1, \ldots, X_n | \theta) f(X_{n+1}|\theta)}{f(X_1, \ldots, X_n | \theta) f(X_{n+1}|\theta)} \right)^2 \right) \\
			&= E \left( \left( \frac{\left(\frac{\partial}{\partial \theta} f(X_1, \ldots, X_n | \theta)\right) f(X_{n+1}|\theta) + f(X_1, \ldots, X_n | \theta) \frac{\partial}{\partial \theta} f(X_{n+1}|\theta)}{f(X_1, \ldots, X_n | \theta) f(X_{n+1}|\theta)} \right)^2 \right) \\
			&= E \left( \left(
			\frac{\frac{\partial}{\partial \theta} f(X_1, \ldots, X_n | \theta)}{f(X_1, \ldots, X_n | \theta)}
			+
			\frac{\frac{\partial}{\partial \theta} f(X_{n+1}|\theta)}{f(X_{n+1}|\theta)}
			\right)^2 \right) \\
			&= E \left( \left(
			\frac{\partial}{\partial \theta} \log f(X_1, \ldots, X_n | \theta)
			+
			\frac{\partial}{\partial \theta} \log f(X_{n+1}|\theta)
			\right)^2 \right) \\
			&= \operatorname{Var}\left( \frac{\partial}{\partial \theta} \log f(X_1, \ldots, X_n | \theta)
			+
			\frac{\partial}{\partial \theta} \log f(X_{n+1}|\theta) \right) \\
			&= \operatorname{Var}\left( \frac{\partial}{\partial \theta} \log f(X_1, \ldots, X_n | \theta)\right)
			+
			\operatorname{Var}\left(\frac{\partial}{\partial \theta} \log f(X_{n+1}|\theta) \right) \\
			&= nI_1(\theta) + I_1(\theta) = (n+1) I_1(\theta) \qedhere
%			I_n(\theta) &= \int_\R \ldots \int_\R \left( \frac{\frac{\partial}{\partial \theta} f(x_1, \ldots, x_n|\theta)}{f(x_1, \ldots, x_n|\theta)}\right)^2 f(x_1, \ldots, x_n|\theta) \, dx_1 \ldots dx_n \\
%			&= \int_\R \ldots \int_\R \frac{\left(\frac{\partial}{\partial \theta} f(x_1, \ldots, x_n|\theta)\right)^2}{f(x_1, \ldots, x_n|\theta)} \, dx_1 \ldots dx_n \\
%			&= \int_\R \ldots \int_\R \frac{\left(\frac{\partial}{\partial \theta} \prod_{i=1}^n f(x_i|\theta)\right)^2}{\prod_{i=1}^n f(x_i|\theta)} \, dx_1 \ldots dx_n \\
%			&= \int_\R \ldots \int_\R \frac{\left(\frac{\partial}{\partial \theta} \prod_{i=1}^n f(x_i|\theta)\right)^2}{\prod_{i=1}^n f(x_i|\theta)} \, dx_1 \ldots dx_n
		\end{align*}
	\end{proof}
	\begin{theorem}
		Let $X_1, \ldots, X_n$ be a random sample. Let $\hat \theta$ be an unbiased estimator of $\theta \in \Theta$. Suppose that $\Theta \subseteq \R$ is an open interval, that $I_n > 0$, as well as
		\[ \int_\R \hat \theta(x) \frac{\partial}{\partial \theta} f(x|\theta) \, dx = \frac{\partial}{\partial \theta} \int_\R \hat \theta(x) f(x|\theta) \, dx \]
		Then
		\[ \operatorname{Var}(\hat \theta) \geq I_n(\theta)^{-1} \]
	\end{theorem}
	We instead discuss this theorem through a more general notion, which has the above theorem as a corollary.
	\begin{theorem}
		Let $g : \Theta \to \R$ be a measurable function. Let $T : \R \to \R$ be an unbiased estimator of $g(\theta)$. Suppose that $\Theta \subseteq \R$ is an open interval, that $I(\theta) > 0$, as well as
		\[ \int_\R T(x) \frac{\partial}{\partial \theta} f(x|\theta) \, dx = \frac{\partial}{\partial \theta} \int_\R T(x) f(x|\theta) \, dx \]
		Then
		\[ \operatorname{Var}(\hat \theta) \geq g'(\theta)^2 I(\theta)^{-1} \]
	\end{theorem}
	\begin{proof}
		Unbiasedness of $T$ gives that $ET = g(\theta)$ for given $\theta \in \Theta$.
		%		We have the following known formula for the Fisher information, established in the lectures:
		%		\[ I(\theta) = -E\left( \frac{\partial^2}{\partial \theta^2} \log f(X|\theta) \right) \]
		%		where $f$ is the density of $X$.
		We compute:
		\begin{align*}
			\operatorname{Var}_\theta(T) I(\theta) &= \left ( \int_\R (T(x) - g(\theta))^2 f(x|\theta) \, dx \right ) \left( \int_\R \left (\frac{\partial}{\partial \theta} \log f(x|\theta)\right )^2 f(x|\theta) \, dx \right) \\
			&\geq \left ( \int_\R (T(x) - g(\theta)) \frac{\partial}{\partial \theta} f(x|\theta) \, dx \right)^2 \\
			&= \left ( \frac{\partial}{\partial \theta} \int_\R T(x) f(x|\theta) \, dx \right)^2 \\
			&= \left ( \frac{\partial}{\partial \theta} g(\theta) \right)^2 \\
			&= g'(\theta)^2 \iff \operatorname{Var}_\theta(T) \geq \frac{g'(\theta)^2}{I(\theta)}
		\end{align*}
		due to the fact that for $f, g \in L^2(\R; \R)$ we have that the following constitutes an inner product:
		\[ \inprod{f, g} := \int_\R f(x) g(x) \, dx \]
		and thus the Cauchy-Schwarz inequality holds, i.e.
		\[ |\inprod{f, g}|^2 \leq \inprod{f, f} \inprod{g, g} \qedhere \]
	\end{proof}
	\newpage
	\subsection{Confidence intervals}
	\begin{definition}
		Let $X_1, \ldots, X_n$ be a random sample with parameter $\theta \in \Theta$. Let $L, U$ be statistics on $X_1, \ldots, X_n$. Let $\alpha \in [0, 1]$ We say that the interval (which is kind of a random variable in its own right) $[L(X_1, \ldots, X_n), U(X_1, \ldots, X_n)]$ is a $1 - \alpha$ \textbf{confidence interval} if
		\[ P(L \leq \theta \leq U) \geq 1 - \alpha \]
		which intuitively means that the parameter does not lie in that interval at most $\alpha$ of the time, so $\alpha$ better be small.\\
		If this inequality is strict, we say that the confidence interval is \textbf{exact}.
	\end{definition}
	\begin{remark}
		Typically, one can compute these confidence intervals by first determining the distance between an unbiased estimator and the parameter, which establishes an interval around the parameter given the estimator. For example, let $X_1, \ldots, X_n \sim N(\mu, \sigma^2)$, we consider the estimator $\bar X \sim N(\mu, \sigma^2/n)$. We want to find an exact confidence interval, which we do by finding $c_\alpha$ such that
		\[ P(|\bar X - \mu| \leq c_\alpha) = 1 - \alpha \]
		We compute
		\begin{align*}
			P(|\bar X - \mu| \leq c_\alpha) &= P(\bar X - \mu \leq c_\alpha) - P(\bar X - \mu \leq -c_\alpha) \\
			&= P(\bar X \leq \mu + c_\alpha) - P(\bar X \leq \mu-c_\alpha) \\
			&= 2P(\bar X \leq \mu + c_\alpha) - 1 \\
			&= 2P((\bar X - \mu)/\sqrt{\operatorname{Var}(\bar X)} \leq c_\alpha/\sqrt{\operatorname{Var}(\bar X))} - 1 = 1 - \alpha
		\end{align*}
		We get that
		\[ c_\alpha := \sqrt{\operatorname{Var}(\bar X))} \Phi^{-1}(1 - \alpha / 2) = \frac{\sigma}{\sqrt{n}} \Phi^{-1}(1 - \alpha / 2) \] 
		does the trick. We get that
		\[ P\left(\bar X-\frac{\sigma}{\sqrt{n}}\Phi^{-1}(1 - \alpha / 2) \leq \mu \leq \bar X+\frac{\sigma}{\sqrt{n}}\Phi^{-1}(1 - \alpha / 2)\right) = 1 - \alpha \]
	\end{remark}
	\begin{remark}
		Confidence intervals could involve any statistics, however, unbiased estimators are obviously more convenient to use.
	\end{remark}
	\newpage
	\subsection{Hypothesis testing}
	\begin{definition}
		Let $\Theta$ be any parameter space. Let $(S, \Sigma)$ be a measure space. Let $\mathcal{P}_{\cdot}(\cdot)~:~\Theta~\to~(\Sigma~\to~[0,~1])$ be a family of probability measures, such that for each $\theta \in \Theta$, we have that $(S, \Sigma, P_\theta)$ is a probability space. Let $H_0, H_1 \subseteq \Theta$ be a disjoint partition of $\Theta$, so $H_0 \cap H_1 = \varnothing$ and $H_0 \cup H_1 = \Theta$. Let $\phi : S \to [0, 1]$ be a measurable function. We say that $\phi$ is a \textbf{hypothesis test} of $H_0$ against $H_1$. We say that it is \textbf{non-randomized} if $\phi(S) = \{0, 1\}$, otherwise we say it is \textbf{randomized}.
	\end{definition}
	\begin{definition}
		We say that a hypothesis $H \subseteq \Theta$ is simple if $|H| = 1$.
	\end{definition}
	\begin{definition}
		We say that for two hypotheses $H_0, H_1$ that the test is
		\begin{itemize}
			\item \textbf{one-sided} if both are intervals in $\R$.
			\item \textbf{two-sided} if $H_0$ is a bounded interval and $H_1$ is its complement.
		\end{itemize}
	\end{definition}
	\begin{remark}
		If we have a simple hypothesis $H = \{\theta_0\}$ for a certain $\theta_0 \in \Theta$, we subscript $P, E, \operatorname{Var}$, etc, with its index ($0$ or $1$), to refer to $P_{\theta_0}$ and the induced $E_{\theta_0}$, etc.
	\end{remark}
	\begin{remark}
		We usually call $H_0$ the \textbf{null hypothesis} and $H_1$ the \textbf{alternative hypothesis}. Of course their roles could be interchanged in the above definitions, however, we usually pick $H_0$ to be the case where there is no observed effect, and $H_1$ where there is (something you are trying to prove). The rest of the subject concerns itself with finding, given a certain confidence $1 - \alpha$, whether $H_0$ is probably false, or nothing can be said.
	\end{remark}
	\begin{example}
		Let $\Theta := \R$ and consider a random variable $X \sim N(\mu \in \Theta, 81)$. Let $H_0 := \{\mu_0\}$, $H_1 := \R \setminus H_0$, which we usually write as $H_0 : \mu = \mu_0$ and $H_1 : \mu \neq \mu_0$. If $H_0$ is true we have that under confidence level $1 - \alpha$, $\hat \mu$ lies in
		\[ C_n := \left[ \mu_0 - \Phi^{-1}(1 - \alpha/2)\frac{9}{\sqrt{n}}, \mu_0 + \Phi^{-1}(1 - \alpha/2)\frac{9}{\sqrt{n}} \right] \]
		that is, $P_{H_0}(\hat \mu \in C_n)$.
		We have the logical implication $H_0 \text{ true} \implies P(\hat \mu \in C_n) = 1 - \alpha$, which kind of also means that $\hat \mu \notin C_n \implies P(H_0) = \alpha$. If indeed you have that $\hat \mu \notin C_n$ we say that we \textbf{reject the null hypothesis} with significance level $\alpha$, since it is very unlikely (for small $\alpha$) that $H_0$ is true.
	\end{example}
	\begin{definition}
		Let $H_0, H_1$ be disjoint hypotheses and let $\phi$ be a non-randomized test on it. Then
		\begin{itemize}
			\item We define the \textbf{power} $\beta(\theta)$ for a parameter $\theta \in \Theta$ by
			\[ \beta(\theta) := P_\theta(\phi = 1) \]
			\item We say that the \textbf{significance level} $\alpha$ is
			\[ \alpha := \sup_{\theta \in H_0} P_\theta(\phi = 1) = \sup_{\theta \in H_0} E_\theta(\phi = 1) \]
			i.e. the least upper bound of the power.
			\item The \textbf{power of the test} for $\theta \in H_1$ is $\beta(\theta)$, and the \textbf{probability of a type I error} (which is when $H_0$ is true yet rejected) is $\beta(\theta)$ for $\theta \in H_0$.
			\item For $\theta \in H_1$, $P_\theta(\phi = 0)$ is the \textbf{probability of a type II error} (which is when $H_0$ is false however it is not rejected)
			\item The \textbf{rejection region} is the set
			\[ \operatorname{RR}_\phi := \{ s \in S \mid \phi(s) = 1 \} \]
		\end{itemize}
	\end{definition}
	\begin{definition}
		Let $T$ be a statistic. Let $P$ be a valence 1 predicate on $\R$. Let $\operatorname{RR} := \{ s \in S \mid P(T(s)) \}$. We also write $T \in \operatorname{RR}$ if $P(T)$ for fixed $s$. We say that $T$ is a \textbf{test statistic} with test $\phi$ such that $\operatorname{RR}_\phi = \operatorname{RR}$.
	\end{definition}
	\begin{definition}
		Let $T$ be a test statistic. Suppose its predicate $P$ is given by $P(T) \iff T \geq c_\alpha$ for critical value $c_\alpha$ at significance level $1 - \alpha$. The $p$\textbf{-value} is defined to be the $\alpha$ such that $T = c_\alpha$ for given (observed) $T$.
	\end{definition}
	\begin{remark}
		This kind of means that the $p$-value is the smallest $\alpha$ for which we would reject the null hypothesis, which means that we reject the null hypothesis at significance levels $1 - \alpha$ if $\alpha \geq p$.
	\end{remark}
	\begin{remark}
		Standard procedure:
		\begin{itemize}
			\item Determine statistical model
			\item Determine hypotheses
			\item Find estimators
			\item Verify conditions
			\item Find distribution of test statistic under null hypothesis
			\item Find critical values / p value for test statistic
			\item Determine significance and report
		\end{itemize}
	\end{remark}
	\begin{example}
		Another example for the $p$-value: let $X_1, \ldots, X_n \sim N(\mu, \sigma^2)$ with known $\sigma$ and test $H_0 : \mu = \mu_0$ vs $H_1 : \mu \neq \mu_0$ for given $\mu_0 \in \R$. Let $\alpha$ be a given significance level. We first construct a confidence interval for $\mu$. We have the estimator $\hat \mu \sim N(\mu, \sigma^2/n)$, so
		\[ Z := \frac{\hat \mu - \mu}{\sigma/\sqrt{n}} \sim N(0, 1) \]
		therefore
		\[ P(|Z| \leq c_\alpha) = 1 - \alpha \iff c_\alpha = \Phi^{-1}(1 - \alpha/2) \]
		So we reject the null hypothesis if $|Z| \geq \Phi^{-1}(1 - \alpha/2)$. $p$-value is thus
		\[ \alpha := 2(1 - \Phi(|Z|)) \]
	\end{example}
	\begin{remark}
		What if $\sigma$ is unknown? Well, we estimate of course. We replace $\sigma$ by $\hat \sigma$. But what distribution do we get?
	\end{remark}
	\begin{definition}
		Let $Z$ be standard normal and let $Y \sim \chi^2_n$ (chi squared with degrees of freedom $n$). We say that
		\[ T_n := \frac{Z}{\sqrt{Y/n}} \]
		has the \textbf{Student's t-distribution} with degrees of freedom $n$.
	\end{definition}
	\begin{remark}
		Why this definition? Well, we are of course interested in the distribution of
		\[ \frac{\hat \mu - \mu}{\sqrt{\hat \sigma^2}/\sqrt{n}} \]
		where
		\[ \hat \sigma^2 := \frac{1}{n-1} \sum_{i=1}^n \left(X_i - \hat \mu\right)^2 \]
		Clearly, $Z := \frac{\hat \mu - \mu}{\sigma/\sqrt{n}}$ is standard normal, and it can be shown that $Y := \hat \sigma^2/\sigma^2$ is $\chi^2_{n-1}$, so indeed we get a Student $t$-distribution. So we can use the CDF of the $t$-distribution, call it $F$, to find $p$-values, compute $p$ for the data in question, and 
	\end{remark}
	\newpage
	\subsection{Neyman-Pearson}
	\begin{definition}
		Let $\alpha \in [0, 1]$ be the significance level. Let $\phi, \psi$ be two statistical tests on the same hypotheses. We define the partial ordering $\leq$ on tests at significance level $\alpha$ by
		\[ \phi \leq \psi \defequiv \alpha_\phi \leq \alpha \land \alpha_\psi \leq \alpha \land \forall \theta \in \Theta_1 : \beta_\phi(\theta) \leq \beta_\psi(\theta) \]
		Alternatively, we say that $\phi$ is \textbf{less powerful} than $\psi$, and vice-versa for \textbf{most powerful}.\\
		A test $\phi$ is said to be \textbf{most powerful} (at a certain level) if for all tests $\psi$ we have $\psi \leq \phi$.
	\end{definition}
	\begin{remark}
		A most powerful test might not exist and might not be unique. The following theorem gives a sufficient condition for when a unique most powerful test exists.
	\end{remark}
	\begin{theorem}
		Let $X$ be a $n$-dimensional random vector of i.i.d. random variables in a probability space $(S, \Sigma, P)$. Let $H_0, H_1$ be two simple hypotheses. Let $L_0, L_1 : S \to [0, \infty)$ be the corresponding densities of $X$. Suppose that $\operatorname{supp}(L_1) \subseteq \operatorname{supp}(L_0)$. Then for all $\alpha \in (0, 1)$ there exist $c \geq 0$ and $y \in [0, 1]$ such that the following test $\phi : X \to [0, 1]$ is most powerful
		\[ \phi(x) := \ind{L_1(x) > cL_0(x)} + \gamma \ind{L_1(x) = cL_0(x)} \]
		Additionally, the test has level $\alpha$, i.e.
		\[ E_0[\phi(X)] = \int_S \phi(x) L_0(x) \, dx = \alpha \]
		and is unique up to choice of $\gamma$.
	\end{theorem}
	\begin{proof}
		Let $\psi$ be another test function. Suppose that it has the same or better level, that is, $E_0[\psi] \leq E_0[\phi]$. This is a fair assumption to make since a test with worse level is never considered to be MP. Let $c \geq 0$ be arbitrary. By cases it can be easily shown that for any $x \in S$
		\[ h(x) := (\phi(x) - \psi(x))(cL_0(x) - L_1(x)) \leq 0 \]
		Namely:
		\begin{itemize}
			\item If $L_1(x) - cL_0(x) > 0$, $\phi(x) - \psi(x) = 1 - \psi(x) \geq 0$ so $h(x) \leq 0$
			\item If $L_1(x) - cL_0(x) < 0$, $\phi(x) - \psi(x) = -\psi(x) \leq 0$ so $h(x) \leq 0$
			\item Otherwise it is trivial
		\end{itemize}
		Integrating yields:
		\[ \int_\R h(x) \, dx = c(E_0(\phi) - E_0(\psi)) - (E_1[\phi] - E_1[\psi]) \leq 0 \]
		so we get the chain
		\[ 0 \leq c(E_0(\phi) - E_0(\psi)) \leq E_1[\phi] - E_1[\psi] \]
		and hence $E_1[\psi] \leq E_1[\phi]$, which proves $\phi$ being MP. Existence of $c$ for given $\alpha$ and a witness $\gamma$ use technical analysis arguments, that I won't bother with.
	\end{proof}
	\begin{remark}
		If the measure of $\{L_1(x) = cL_0(x)\}$ is zero, i.e. it has probability zero, i.e. they do not agree a.e., $\gamma$ can be picked arbitrarily in $[0, 1]$, so we can either choose to reject on the boundary ($\gamma = 1$) and other stuff, or introduce randomization.
	\end{remark}
	
	\newpage
	\section{Regression}
	Let $X_1, \ldots, X_n$ be random variables (`the input') and $Y$ as well (`the output'). We want to find out whether there exists some noise variable $\varepsilon$ and a measurable function $f : \R^n \to \R$ such that $Y = f(X_1, \ldots, X_n) + \varepsilon$, i.e. whether such a relation exists. This turns out to be quite the impossible question, so we make simplifications. We restrict to linear $f$. Not linear in the vector space sense, but in the sense that there exists a vector space linear $g : \R^n \to \R$ such that $f(x_1, \ldots, x_n) = g(x_1, \ldots, x_n) + b$ for $b \in \R$, i.e. this is an \textbf{affine} $f$. We essentially want to solve
	\[ Y = \beta_0 + \beta_1 X_1 + \ldots + \beta_n X_n + \varepsilon \]
	with $\beta_0$ the \textbf{intercept}, \textbf{regression coefficients} $\beta_1, \ldots, \beta_n$ and \textbf{noise} $\varepsilon$. This is called a \textbf{linear model}.
	\subsection{Preliminaries}
	\begin{definition}
		Let $X$ and $Y$ be random variables. We say that the quantity $\rho_{X, Y}$ is the \textbf{correlation coefficient} of $X$ and $Y$, where
		\[ \rho_{X, Y} := \frac{\operatorname{Cov}(X, Y)}{\sqrt{\operatorname{Var}(X)\operatorname{Var}(Y)}} \]
		where
		\[ \operatorname{Cov}(X, Y) := E\left[(X - EX)(Y - EY)\right] \]
		It can be shown that $|\rho_{X, Y}| \leq 1$, and quite trivially that for all $a, b \in \R \setminus \{0\}$
		\[ \rho(aX, bY) = \rho(X, Y) \]
		Also clearly when $X$ and $Y$ are independent, then $\rho_{X, Y} = 0$.
	\end{definition}
	\begin{definition}
		Let $X, Y$ be random variables. We say that $X$ and $Y$ are \textbf{uncorrelated} if $\rho_{X, Y} = 0$.
	\end{definition}
	\begin{remark}
		Both the covariance and the correlation coefficient can be computed from the sample, simply by using estimators, e.g.
		\[ \operatorname{\hat{Cov}}(X, Y) := \frac{1}{n} \sum_{i=1}^n (X_i - \bar X)(Y_i - \bar Y) \]
		\[ \hat{\rho}(X, Y) := \frac{\sum_{i=1}^n (X_i - \bar X)(Y_i - \bar Y)}{\sqrt{\sum_{i=1}^n (X_i - \bar X)^2 \sum_{i=1}^n(Y_i - \bar Y)^2}} \]
	\end{remark}
	\newpage
	\subsection{Covariance matrix -- Computing linear model from theory}
	\begin{definition}
		Let $V$ be a finite-dimensional inner product vector space. Let $X$ be a random variable on $V$. Let $\mu := EX$ be its mean. We define the map $\Sigma_X : V \to V$ by first establishing that for every vector $v \in V$ there is a unique dual $v^* \in V^*$ such that for all $x \in V$, $v^*(x) = \inprod{x, v}$, i.e. $v^* := \inprod{\cdot, v}$. We define then for $v \in V$
		\[ \Sigma_X(v) := E[(X - \mu)(X - \mu)^*(v)] \]
		This is usually written more compactly as
		\[ \Sigma_X := E[(X - \mu)(X - \mu)^*] \]
		or, when working in an orthonormal basis exclusively
		\[ \Sigma_X := E[(X - \mu)(X - \mu)^T] \]
		Indeed, if we take components in an orthonormal basis $e_1, \ldots, e_n$ we get
		\begin{align*}
			(\Sigma_X)^i_j &= \varepsilon^i(E[(X - \mu)(X - \mu)^*(e_j)]) \\
			&= \varepsilon^i(E[(X - \mu)(X - \mu)^j]) \\
			&= E[(X - \mu)^i(X - \mu)^j] = \operatorname{Cov}(X_i, X_j)
		\end{align*}
	\end{definition}
	\begin{lemma}
		Very easily, we have that
		\begin{itemize}
			\item In an orthonormal basis, $\Sigma_X$ is symmetric
			\item $\Sigma_X$ is positive semi-definite, that is, for all $v \in V$, we have $v^*(\Sigma_X(v)) \geq 0$.
		\end{itemize}
	\end{lemma}
%	\begin{theorem}
%		Let $X_1, \ldots, X_n$ and $Y$ have a linear model with noise $\varepsilon$. Suppose that $X = (X_1, \ldots, X_n)$ has invertible covariance matrix and the noise $\varepsilon$ is uncorrelated with $X$. Then
%		\[ \Sigma_X \beta = \operatorname{Cov}(X, Y) := (\operatorname{Cov}(X_1, Y), \ldots, \operatorname{Cov}(X_n, Y)) \]
%		
%	\end{theorem}
	\begin{theorem}
		Let $X_1, \ldots, X_n$ and $Y$ have a linear model with noise $\varepsilon$. Suppose the noise $\varepsilon$ is uncorrelated with $X$. Then
		\[ \Sigma_X \beta = \operatorname{Cov}(X, Y) := (\operatorname{Cov}(X_1, Y), \ldots, \operatorname{Cov}(X_n, Y)) \]
		where $\beta := (\beta_1, \ldots, \beta_n)$.
	\end{theorem}
	\begin{proof}
%		We first compute:
%		\[ EY = \beta_0 + E\inprod{\beta, X} + E\varepsilon \]
		By components:
		\begin{align*}
			(\Sigma_X \beta)^i &= E[(X - \mu)^i \inprod{\beta, X - \mu}] \\
			&= E\left[(X - \mu)^i \left(\beta_1(X_1 - EX_1) + \ldots + \beta_n(X_n - EX_n)\right)\right] \\
			&= \operatorname{Cov}(X_i, \beta_1 X_1 + \ldots + \beta_n X_n) \\
			&= \operatorname{Cov}(X_i, Y - \beta_0 - \varepsilon) \\
			&= \operatorname{Cov}(X_i, Y) \qedhere
%			&= E\left[(X - \mu)^i (Y - \beta_0 - \varepsilon - \inprod{\beta, \mu})\right] \\
%			&= E\left[(X - \mu)^i (Y - \beta_0 - E\inprod{\beta, X} - E\varepsilon)\right] \\
%			&= E\left[(X - \mu)^i (Y - EY)\right] \\
%			&= \operatorname{Cov}(X_i, Y)
		\end{align*}
	\end{proof}
	\begin{corollary}
		Suppose that $X = (X_1, \ldots, X_n)$ has invertible covariance matrix and the noise $\varepsilon$ is uncorrelated with $X$, then
		\[ \beta = \Sigma_X^{-1} \operatorname{Cov}(X, Y) \]
		and
		\[ \beta_0 = EY - \inprod{EX, \beta} - E\varepsilon \]
		Typically we assume that $E\varepsilon = 0$ so $\beta_0 = EY - \inprod{EX, \beta}$.
	\end{corollary}
	\begin{corollary}
		If $\Sigma_X$ is symmetric positive definite, it is invertible, so we also have the above.
	\end{corollary}
	\begin{proof}
		Classic linear algebra exercise. One can show that symmetric positive definite matrices have strictly positive eigenvalues, so it must be invertible; otherwise, $0$ is an eigenvalue.
	\end{proof}
	\begin{corollary}
		To check if the covariance matrix is invertible, it suffices to check whether $0$ is an eigenvalue.
	\end{corollary}
	\newpage
	\subsection{MSE estimation -- Computing linear model from data}
	\begin{definition}
		Let $(x_{i, 1}, \ldots, x_{i, n}, y_i)$ for $i = 1, \ldots, m$ be $m$ samples / data of a linear model:
		\[ y_i = \beta^j x_{i, j} + \beta_0 + \varepsilon^i \]
		where we used the Einstein convention. We can now compute the error $\varepsilon^i$, and we define the \textbf{mean squared residual error} as
		\[ \operatorname{MSE}(\beta) := \frac{1}{n} \sum_{i=1}^n (\varepsilon^i)^2 \]
	\end{definition}
	\begin{remark}
		If we assume that $E\varepsilon = 0$, it is reasonable to expect that a good model minimizes at least the second moment of $\varepsilon$, which means that $\operatorname{MSE}(\beta)$ is minimized. But, we already know how to solve this exact problem, from linear algebra. Observe that, using the standard inner product and its induced norm, we have that
		\[ \operatorname{MSE}(\beta) = \frac{1}{n}\norm{y - X\beta}^2 \]
		where
		\[ X = \bmat{
			1 & x_{1, 1} & \ldots & x_{1, n} \\
			\vdots & & \vdots & \\
			1 & x_{m, 1} & \ldots & x_{m, n}
		} \]
		If $\im X = \R^m$, we have that $X^TX$ is invertible, so a nasty computation shows that
		\[ X^T(X \beta - y) = 0 \]
		minimizes $\operatorname{MSE}(\beta)$, so we get
		\[ \beta = (X^TX)^{-1}(X^Ty) \]
	\end{remark}
	\subsection{Application in $\R^2$}
	We consider the linear model
	\[ Y = \beta_0 + \beta_1 X + \varepsilon \]
	Now we get that
	\[ \beta_1 = \frac{\operatorname{Cov}(X, Y)}{\operatorname{Cov}(X, X)}, \qquad \beta_0 = EY - \beta_1 EX \]
	We can use this on data with the sample covariance $\operatorname{\hat{Cov}}$ and sample standard deviation $\hat \sigma$ and of course the sample means. The MSE is also easy to compute:
	\[ \operatorname{MSE}(\beta) = \operatorname{Var}(\varepsilon) = \operatorname{Var}(Y - \beta_1 X) = \operatorname{Var}(Y) + \beta_1^2\operatorname{Var}(X) - 2\beta_1 \operatorname{Cov}(X, Y) \]
	This is minimized whenever
	\[ \beta_1 = \frac{\operatorname{Cov}(X, Y)}{\operatorname{Var}(X)} \]
	which, surprise surprise, is the same as the theoretical `optimal'. So, under the assumption that data follows a linear model, the MSE approximation is pretty dang good.
	
	\newpage
	\subsection{Gaussian noise}
	We make the same assumptions as before, as well as $\varepsilon \sim N(0, \sigma^2)$ for known $\sigma$. We already know how to deal with this, with the least-squares estimator. But we can also find an MLE for $\beta_0, \beta_1$. Assume we have random samples $Y_1, \ldots, Y_n$ according to
	\[ Y_i = \beta_0 + \beta_1 x_i + \varepsilon_i \]
	with fixed $x_1, \ldots, x_n$. We clearly have that $Y_i \sim N(\beta_0 + \beta_1 x_i, \sigma^2)$. So
	We compute the likelihood
	\begin{align*}
		f_{Y_1, \ldots, Y_n}(y_1, \ldots, y_n|\beta_0, \beta_1, \sigma^2) &= \prod_{i=1}^n \frac{1}{\sqrt{2 \pi \sigma^2}} \exp\left( - \frac{(y_i - (\beta_0 + \beta_1 x_i))^2}{2 \sigma} \right) \\
		&= (2 \pi \sigma^2)^{-n/2} \exp\left( - \sum_{i=1}^n \frac{(y_i - (\beta_0 + \beta_1 x_i))^2}{2 \sigma} \right)
	\end{align*}
	The log-likelihood is
	\[ \log f_{Y_1, \ldots, Y_n}(y_1, \ldots, y_n|\beta_0, \beta_1, \sigma^2) = \ldots - \frac{1}{2\sigma} \sum_{i=1}^n (y_i - \beta_0 - \beta_1 x_i)^2 \]
	which is clearly maximized whenever $\sum_{i=1}^n (y_i - \beta_0 - \beta_1 x_i)^2$ is minimized, which recovers the least-squares problem and indeed we get the same solution as before. It is a bit more work to prove that then also the $\hat \sigma^2$ etc estimators maximize the likelihood, but not that difficult to see. We define
	\[ \hat Y_i := \hat \beta_0 + \hat \beta_1 x_i \]
	Now when we have this Gaussian assumption we say that the measured noise values are the residuals, denoted by $R_i := \hat Y_i - \hat \beta_0 - \hat \beta_1 x_i$, which should in theory be $N(0, \sigma^2)$ distributed.
	\begin{lemma}
		We also have that in this model
		\[ \sum_{i=1}^n (x_i - \bar x)(Y_i - \bar Y) = \sum_{i=1}^n (x_i - \bar x)Y_i \]
	\end{lemma}
	\begin{proof}
		\begin{align*}
			\sum_{i=1}^n (x_i - \bar x)(Y_i - \bar Y)
			&= \sum_{i=1}^n (x_i - \bar x)Y_i + \bar Y\underbrace{\sum_{i=1}^n (x_i - \bar x)}_{\sum_i X_i - n\bar x = 0} \\
			&= \sum_{i=1}^n (x_i - \bar x)Y_i \qedhere
		\end{align*}
	\end{proof}
	\begin{lemma}
		We have that (where $\sigma^2 = \operatorname{Var}(\varepsilon_i)$)
		\[ \hat \beta_1 \sim N\left(\beta_1, \frac{\sigma^2}{\sum_i (x_i - \bar x)^2}\right) \]
	\end{lemma}
	\begin{proof}
		Trivially we have $E\hat \beta_1 = \beta_1$ and
		\[ \operatorname{Var}(\hat \beta_1) = \frac{\operatorname{Var}(Y_1) \sum_{i=1}^n (x_i - \bar x)^2}{\left(\sum_{i=1}^n (x_i - \bar x)^2\right)^2} = \frac{\sigma^2}{\sum_i (x_i - \bar x)^2} \qedhere \]
	\end{proof}
	
	\newpage
	\subsection{Unknown variance}
	\begin{lemma}
		The MLE of the Gaussian model from before is (partially repeated)
		\begin{align*}
			\hat \beta_1 &= \frac{\sum_{i=1}^n (x_i - \bar x)Y_i}{\sum_{i=1}^n (x_i - \bar x)^2} = \sqrt{\frac{\hat \sigma^2_Y}{\hat \sigma^2_X}} \hat \rho(X, Y) \\
			\hat \beta_0 &= \bar Y - \hat \beta_1 \bar x
		\end{align*}
		which are unbiased.
	\end{lemma}
	So the above is a good way of computing our desired result, but there is one issue, if we want to study residuals we better know $\sigma$.
	\begin{lemma}
		The MLE for $\sigma$ is
		\[ \hat \sigma := \frac{1}{n} \sum_{i=1}^n (Y_i - \hat \beta_0 - x_i \hat \beta_1)^2 \]
		An unbiased, but with less likelihood, estimator, is
		\[ \tilde \sigma := \frac{1}{n - 2} \sum_{i=1}^n (Y_i - \hat \beta_0 - x_i \hat \beta_1)^2 \]
	\end{lemma}
	\subsection{Hypothesis test}
	We are interested in knowing whenever $\beta_1 \neq 0$, since that tells you that there is at least some correlation, i.e. the correlation coefficient is nonzero. We can do this with hypothesis testing.\\
	We test $H_0 : \beta_1 = 0$ versus $H_1 : \beta_1 \neq 0$. Still under the Gaussian noise assumption, we have that
	\[ \hat \beta_1 \sim N\left(\beta_1, \frac{\sigma^2}{\sum_i (x_i - \bar x)^2} \right) \]
	We also define
	\[ S_{xx} := \sum_{i=1}^n (x_i - \bar x)^2 \]
	such that we have $Z := \frac{\hat \beta_1 - \beta_1}{\sigma/\sqrt{S_{xx}}} \sim N(0, 1)$. Let $\alpha$ be the significance level, then we have a $1 - \alpha$ confidence interval for $Z$, namely $P\left(|Z| \leq \frac{\sigma}{\sqrt{S_{xx}}} \Phi^{-1}(1 - \alpha/2) \right) = 1 - \alpha$, so we reject if \ldots.
	If $\sigma$ is unknown we use the $\hat \sigma^2$ estimator which provably makes $Z$ $t$ distributed with $n-2$ dof.
	
	\newpage
	\subsection{Error}
	If $Y_i$ follows a linear model and $\hat Y_i$ is an estimator / prediction, we can quantify exactly how far off our prediction is given data. We define
	\begin{align*}
		SS_{error} &:= \sum_{i=1}^n (Y_i - \hat Y_i)^2 \\
		SS_{regr} &:= \sum_{i=1}^n (\hat Y_i - \bar Y)^2 \\
		SS_{total} &:= \sum_{i=1}^n (Y_i - \bar Y_i)^2 = SS_{error} + SS_{regr} \\
	\end{align*}
	The first quantifies how far off the prediction is for these specific measurements, the second quantifies how far off the prediction is of the mean, and the last measures to what extent the data even follows a linear pattern, which turns out to be independent of choice of regression. We also define
	\[ R^2 := \frac{SS_{regr}}{SS_{total}} = 1 - \frac{SS_{error}}{SS_{total}} = \hat \rho(X, Y)^2 \]
	So we can already find the correlation coefficient simply by testing how well the regression fits, neat.
	
	\newpage
	\section{$\chi^2$ tests}
	\subsection{Goodness of fit}
	Let $X_1, \ldots, X_n$ be categorical variables with range $\{1, \ldots, m\}$ for $m$ categories. We want to test whether these categories were distributed `fairly' from data, such that category $i$ has probability $p_i$ and we then also have $\sum_i p_i = 1$. This is our null hypothesis. For every $X_j$ we clearly have $P(X_j = i) = p_i$ by setup. We can test this hypothesis by first making a histogram, i.e. defining the tuple $N = (N_1, \ldots, N_m)$ where $N_i$ is defined to be the amount of times any observation equals category $i$. We have that
	\begin{align*}
		P(N = (n_1, \ldots, n_m)) &= P(N_1 = n_1, \ldots, N_m = n_m) \\
		&= \ind{\sum_i n_i = n} \frac{n!}{\prod_i n_i!} \prod_{i=1}^m p_i^{n_i}
	\end{align*}
	This is a generalization of the binomial distribution, called the \textbf{multinomial distribution}, denoted by $N \sim \operatorname{NM}(n, p_1, \ldots, p_m)$.
	Here it does not really make sense to talk about the variance of $N$ but we can find its covariance matrix. We have also that for large $n$, $N$ follows a multivariate normal distribution. Using a bunch of talk that we do not really need to know about we get that
	\[ \sum_{i=1}^m \frac{(N_i - np_i)^2}{np_i} \sim \chi^2_{m-1} \]
	approximately! And this guy should be small under the null hypothesis, so in other words, if the above quantity is bigger than critical values of $\chi^2$, we know that the true proportions likely differ.
	\newpage
	\part{Mathematical Statistics II}
	\section{Random variables on vector spaces}
	Let $V$ be a finite-dimensional $\R$-vector space that is also integrable. A \textbf{random variable} $X$ on $V$ and a probability space $(S, \Sigma, P)$ is a map $X : S \to V$. If $X$ has a density $f : V \to \R$, we define the \textbf{expectation} of $X$ through
	\[ EX := \int_V x f(x) \, dx \]
	which is well-defined due to all densities being equal up to $L^2$ quotient.
	\begin{lemma}
		We have that (in any basis) $(EX)^i = E(X^i)$, where $X^i : S \to \R$ is the $i$-th component of $X$ (pointwisely).
	\end{lemma}
	\begin{lemma}
		For any two random variables $X, Y$ on the same measure and vector spaces we have $E(X + Y) = EX + EY$. That is, $E$ is a linear operator from the space of random variables to $V$. The space of random variables is also a vector space (actually a module), where the field (ring) is the random variables $S \to \R$.
	\end{lemma}
	\begin{proof}
		Go into a basis and establish the equality componentwise through elementary results of probability theory.
	\end{proof}
	\begin{lemma}
		Let $A : V \to W$ be a linear map on vector spaces $V, W$. Then it holds that for any random variable $X$ we have $E[AX] = AE[X]$.
	\end{lemma}
	\begin{proof}
		Trivial due to linearity. We can pick a basis $e_1, \ldots, e_n$ for $V$ for example to have
		\[ E[AX] = A(e_i) E[X^i] = A(e_i) (EX)^i = A(EX) \]
	\end{proof}
	\begin{definition}
		The \textbf{covariance operator} $\operatorname{Cov}(X)$ for a random variable $X$ is given in a basis through
		\[ \operatorname{Cov}(X)^i_j := \operatorname{Cov}(X^i, X^j) \]
	\end{definition}
	\begin{lemma}
		On an inner product space with an orthonormal basis, we also have that
		\[ \operatorname{Cov}(X) = E[\inprod{X - EX, X - EX}] \]
		where the inner product is understood componentwise I suppose.
	\end{lemma}
	\begin{lemma}
		Let $A : V \to W$ be a linear map and let $X$ be a random variable on $V$. We have that
		\[ \operatorname{Cov}(AX) = A\operatorname{Cov}(X)A^* \]
		where $A^*$ is the Hermitian adjoint. Remember that this is the unique linear map such that for any two vectors $x, y \in V$ we have $\inprod{A(x), y} = \inprod{x, A^*(y)}$.
	\end{lemma}
	\begin{proof}
		We pick an orthonormal basis and we work it out. Remember that on an orthonormal basis we have the conjugate transpose blabla. Observe that
		\begin{align*}
			\operatorname{Cov}(AX)^i_j &= \operatorname{Cov}(AX^i, AX^j) \\
			&= \operatorname{Cov}(A^i_k X^k, A^j_k X^k) \\
			&= A^i_k A^j_k \operatorname{Cov}(X^k, X^k) \\
			&= A^i_l (A^*)^k_j \operatorname{Cov}(X^l, X^k) \\
			&= A^i_l \operatorname{Cov}(X)^l_k (A^*)^k_j \\
			&= (A\operatorname{Cov}(X))^i_k (A^*)^k_j \\
			&= (A\operatorname{Cov}(X)A^*)^i_j
		\end{align*}
	\end{proof}
	\begin{definition}
		A linear map $A : V \to W$ on inner product spaces is said to be \textbf{positive semi-definite} is for all $v \in V$ we have $\inprod{v, A(v)} \geq 0$.
	\end{definition}
	\begin{lemma}
		The covariance operator is positive semi-definite.
	\end{lemma}
	\begin{definition}
		A random variable $X$ is said to be \textbf{multivariate normally distributed} with mean $\mu \in V$ and invertible covariance matrix $\Sigma$ if	its density $f_X : V \to \R$ satisfies
		\[ f_X(x) = \frac{1}{\sqrt{\det(2 \pi \Sigma)}} \exp\left( -\frac{1}2 \inprod{x - \mu, \Sigma^{-1}(x - \mu)} \right) \]
		and we write $X \sim N(\mu, \Sigma)$.
	\end{definition}
	\begin{lemma}
		If $X \sim N(\mu, \Sigma)$ then $EX = \mu$ and $\operatorname{Var}(X) = \Sigma$.
	\end{lemma}
	\begin{remark}
		You don't necessarily have that if $X$ and $Y$ are normal, then $(X, Y)$ is multivariate normal. Example, let $X \sim N(0, 1)$. Let $\varepsilon$ be independent noise, any r.v, i.e. $\operatorname{Cov}(X, \varepsilon) = 0$, such that also $P(\varepsilon = 1) = P(\varepsilon = -1) = 1/2$. Let $Y = \varepsilon X$. We can show that $X$ and $Y$ are distributionally equal, i.e. $P(Y \leq y) = P(X \leq y)$. Then the support of the standard normal is not all of $V$, namely, only the subset where $Y = X$ or $Y = -X$.
	\end{remark}
	\begin{lemma}
		If $X$ and $Y$ are normally distributed and independent, then together they are bivariate normal.
	\end{lemma}
	\begin{proof}
		The joint density can be computed. Let $X \sim N(\mu, \sigma^2)$ and $Y \sim N(\tilde \mu, \tilde \sigma^2)$. Let $\mu' := (\mu, \tilde \mu)$ and
		\[ \Sigma := \begin{bmatrix}
			\sigma^2 & 0\\
			0 & \tilde \sigma^2
		\end{bmatrix}, \qquad \Sigma^{-1} = \begin{bmatrix}
			1/\sigma^2 & 0\\
			0 & 1/\tilde \sigma^2
		\end{bmatrix} \]
		Then we have $(X, Y) \sim N(\mu', \Sigma)$. We can see this by
		\begin{align*}
			f_{X, Y}(x, y) &= f_X(x) f_Y(y) \\
			&= \frac{1}{\sqrt{2 \pi \sigma^2} \sqrt{2 \pi \tilde \sigma^2}} \exp\left(- \frac{1}{2} \left( \left(\frac{x - \mu}{\sigma}\right)^2 + \left(\frac{y - \tilde \mu}{\tilde \sigma}\right)^2 \right) \right) \\
			&= \frac{1}{\sqrt{\det(2 \pi \Sigma)}} \exp\left( -\frac{1}2 \inprod{(x, y) - \mu', \Sigma^{-1}((x, y) - \mu')} \right) \qedhere
		\end{align*}
	\end{proof}
	\begin{corollary}
		Let $X_1, \ldots, X_n \sim N(\mu, \sigma^2)$ be i.i.d. Then $X := (X_1, \ldots, X_n)$ is multivariate normal with $\mu' := (\mu, \ldots, \mu)$ and $\Sigma := \sigma^2 \id$ so that $X \sim N(\mu', \Sigma)$.
	\end{corollary}
	\begin{lemma}
		Let $V$ and $W$ be appropriate vector spaces and let $X$ be a random variable on $V$. Let $A : V \to W$ be linear and invertible. Let $Y := AX$. We have that
		\[ f_Y(y) = f_X(A^{-1}(y))|\det(A^{-1})| \]
	\end{lemma}
	\begin{proof}
		The density is the unique function that satisfies that for all measurable $U \subseteq W$ we have
		\[ P(Y \in U) = \int_U f_Y(y) \, dy \]
		So we compute the above:
		\begin{align*}
			P(Y \in U) &= P(AX \in U) \\
			&= P(X \in A^{-1}(U)) \\
			&= \int_{A^{-1}(U)} f_X(x) \, dx \\
			&= \int_U f_x(A^{-1}(y)) |\det(A^{-1})| \, dy \qedhere
		\end{align*}
	\end{proof}
	\begin{lemma}
		Let $V$ and $W$ be appropriate vector spaces and let $X \sim N(\mu, \Sigma)$ on $V$. Let $A : V \to W$ be linear and invertible. We have that
		\[ AX \sim N(A\mu, A \circ \Sigma \circ A^*) \]
		where $A^*$ is the Hermitian adjoint. Recall that this is the unique map that satisfies
		\[ \inprod{A(x), y} = \inprod{x, A^{*}(y)} \]
	\end{lemma}
	\begin{proof}
		We go into components in an orthonormal basis. We then compute the density. Let $Y := AX$.
		\begin{align*}
			f_{Y}(y) &= f_X(A^{-1}(y)) |\det(A^{-1})| \\
			&= \frac{1}{\sqrt{\det(2 \pi \Sigma)}} \exp\left( -\frac{1}2 \inprod{A^{-1}(y) - \mu, \Sigma^{-1}(A^{-1}(y) - \mu)} \right) |\det(A^{-1})| \\
			&= \frac{1}{\sqrt{\det(2 \pi \Sigma)}} \exp\left( -\frac{1}2 \inprod{A^{-1}(y - A\mu), (A^* \circ (A^*)^{-1} \circ \Sigma^{-1} \circ A^{-1})(y - A\mu)} \right) |\det(A)|^{-1} \\
			&= \frac{1}{\sqrt{\det(2 \pi A^2 \Sigma)}} \exp\left( -\frac{1}2 \inprod{y - A\mu, (A \circ \Sigma \circ A^*)^{-1}(y - A\mu)} \right) \\
			&= \frac{1}{\sqrt{\det(2 \pi A^2 \Sigma)}} \exp\left( -\frac{1}2 \inprod{y - A\mu, (A \circ \Sigma \circ A^*)^{-1}(y - A\mu)} \right)
		\end{align*}
		We are done if we can show that $\det(A^2 \circ \Sigma) = \det(A \circ \Sigma \circ A^*)$. Well, we have that $\det(A^*) = \det(A)$ due to previous linear algebra results, and because everything is invertible this is clear.
	\end{proof}
	\begin{lemma}
		Let $X_1, \ldots, X_n \sim N(\mu, \sigma^2)$ be independent. Then
		\[ \frac{(n-1) \hat \sigma^2 }{\sigma^2} = \frac{1}{\sigma^2} \sum_{i=1}^n (X_i - \bar X)^2 \sim \chi^2_{n-1} \]
		where $\bar X := \frac{1}{n} \sum_{i=1}^n X_i$ and $\hat \sigma^2 := \frac{1}{n-1} \sum_{i=1}^n (X_i - \bar X)^2$.
	\end{lemma}
	\begin{proof}{\textbf{TODO}}
		Let $X := (X_1, \ldots, X_n)$.
		\[ (n-1)\hat \sigma^2 = \inprod{X - \bar X, X - \bar X} \qedhere \]
	\end{proof}
	\newpage
	\section{Regression (continued)}
	\subsection{Simple linear regression}
	Let $\beta_0 \in \R$ and $\beta \in \R^n$, $\varepsilon_1, \ldots, \varepsilon_m$ be independent noise with mean $0$, let additionally for $i=1, \ldots, m$ we have $X_i$ be random variables in $\R_n$. We say that $Y := (Y_1, \ldots, Y_m)$ follows a \textbf{linear model} if
	\[ Y_i = \beta_0 + \beta_j X_i^j + \varepsilon_i \]
	Let
	\[ X := \begin{bmatrix}
		1 & X_1^1 & \cdots & X_1^n \\
		\vdots & \vdots & \ddots & \vdots \\
		1 & X_m^1 & \cdots & X_m^n \\
	\end{bmatrix}, \beta := \bmat{\beta_0 \\ \vdots \\ \beta_n} \]
	so that we also have $Y = X\beta + \varepsilon$.
	\begin{lemma}
		If $\varepsilon \sim N(0, \sigma^2 \id)$ then $Y \sim N(X \beta, \sigma^2 \id)$ (where $X_i$ are not assumed to be random variables but concrete data).
	\end{lemma}
	\begin{theorem}
		If $X^* \circ X$ is invertible, then the MLE and the least-squares estimator $\hat \beta$ for $\beta$ coincide and
		\[ \hat \beta = (X^* \circ X)^{-1}(X^*(Y)) \]
		with
		\[ \hat \beta \sim N(\beta, \sigma^2 (X^* \circ X)^{-1}) \]
	\end{theorem}
	\begin{proof}
		The likelihood of $Y$ is
		\[ f(x) = \frac{1}{\sqrt{2 \pi \sigma^2}} \exp\left( -\frac{1}2 \left(\frac{\norm{x - X \beta}}{\sigma}\right)^2 \right) \]
		This is maximized by
		\[ \hat \beta := \operatorname*{argmin}_{\beta} \norm{Y - X \beta} \]
		which is exactly the least-squares estimator given by solving the normal equation. It is unbiased due to
		\[ E\hat \beta = \beta + E(X^* \circ X)^{-1}(X^*(\varepsilon)) = \beta \]
		and normality and the other properties can be derived from the earlier lemma's.
	\end{proof}
	\begin{remark}
		If $X^* X$ is not invertible, we can interpret that as the minimum being attained by infinitely many $\beta$, or alternatively that there are simply multiple linear models possible given the data. Such a model is called \textbf{non-identifiable}.
	\end{remark}
	\newpage
	\subsection{Logistic regression}
	Behold, a majestic smooth interpolator and actually bijection between $(0, 1)$ and $\R$: $f : \R \to (0, 1)$ defined by $x \mapsto f(x) := \left(1 + \exp(-x)\right)^{-1}$. So if we want to do a regression on a categorical variable $X \in \{A, B\}$, we can simply model the probability of category $A$, by taking a linear model (something we understand well) and stuffing it into the logistic function, to get a probability out. I.e. we use the model
	\[ P(X = A) = \left(1 + \exp\left(-\left(\beta_0 + \beta_1 x_1 + \ldots + \beta_n x_n\right)\right)\right)^{-1} \]
	Now let's invert the logistic function, which is by the way invertible as we saw. For $x \in \R$ and $y \in (0, 1)$ we have
	\[ y = (1 + \exp(-x))^{-1} \iff x = -\log(y^{-1} - 1) = \log\left(\frac{y}{1-y}\right) =: \operatorname{logit}(y) \]
	i.e. we define the \textbf{logit} function to be the inverse of the logistic function. We are left with
	\[ \operatorname{logit}(P(X = A)) = \beta_0 + \beta_1 x_1 + \ldots + \beta_n x_n \]
	If we want to now make an estimator, we got to have a statistical model. We still observe data like before and assume that the outputs $Y_i$ take either $0$ or $1$ with probability $p_i$, i.e. Bernoulli distributed.
	Let us compute the likelihood
	\[ L = \prod_{i=1}^n (1-p_i)^{1 - Y_i} p_i^{Y_i} \]
	so that the log-likelihood is
	\[ \log L = \sum_{i=1}^n (1 - Y_i) \log(1-p_i) + Y_i \log(p_i) \]
	where
	\[ p_i = 1/(1 + \exp(\beta_0 + \beta_1 x_{i,1} + \ldots + \beta_n x_{i,n})) \]
	well, good luck maximizing that\textellipsis
	
\end{document}
