{ "nbformat_minor": 2, "nbformat": 4, "cells": [ { "source": [ "$$\n", "\\def\\CC{\\bf C}\n", "\\def\\QQ{\\bf Q}\n", "\\def\\RR{\\bf R}\n", "\\def\\ZZ{\\bf Z}\n", "\\def\\NN{\\bf N}\n", "$$\n", "# Introduction to Markov chains\n", "\n", "In this worksheet we simulate a random process called *Markov chain*. Our source of randomness will be the function `random()` that returns a random floating point number between 0 and 1 uniformly" ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "print random()" ], "outputs": [], "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "print random()" ], "outputs": [], "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "[random() for i in range(10)]" ], "outputs": [], "metadata": {} }, { "source": [ "## Discrete random variables from `random()`\n", "\n", "### Coin tossing\n", "\n", "Write a Sage function `toss()` that returns $0$ or $1$, each with probability $1/2$." ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "### Biased coin toss (a.k.a. Bernouilli measure)\n", "\n", "Write a Sage function `biased_toss(p)` that, given a number $p\\in\n", "[0,1]$ returns $0$ with probability $p$ and $1$ with probability $1-p$." ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "### Finite probabilities\n", "\n", "Let $p = (p_0, p_1, \\ldots, p_{n-1}) \\in [0,1]^n$ be such that $\\sum_{j=0}^{n-1} p_j = 1$. Write a function `multi_toss(p)` so that given a vector $p$ returns an integer $\\{0, 1, \\ldots, n-1\\}$ with probability $p_j$." ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "## Markov chain simulation\n", "\n", "Now we consider a *stochastic matrix*, ie a square matrix $P = (p_{ij})_{0 \\leq i,j < n}$ with non-negative coefficients and such that\n", "\n", "$$\\forall i \\in \\{0, 1, \\ldots, n-1\\}, \\quad \\sum_{j=0}^{n-1} p_{ij} = 1$$\n", "\n", "Here are four examples" ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "P1 = 1/3 * matrix(4, [0,3,0,0,1,0,2,0,0,2,0,1,0,0,3,0])\n", "P2 = 1/12 * matrix(5, [0,4,4,4,0,4,0,0,4,4,6,0,0,6,0,3,3,3,0,3,0,6,0,6,0])\n", "P3 = 1/10 * matrix(5, [6,4,0,0,0,3,7,0,0,0,0,5,0,5,0,0,0,0,2,8,0,0,0,6,4])\n", "P4 = 1/3 * matrix(4, [0,1,1,1,1,0,1,1,1,1,0,1,1,1,1,0])" ], "outputs": [], "metadata": {} }, { "source": [ "Write a function `is_stochastic_matrix(P)` that checks whether the input $P$ is a stochastic matrix and returns `True` or `False` accordingly" ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "Check that the matrices `P1`, `P2`, `P3` and `P4` are stochastic matrices" ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "Let $P = (p_{ij})_{0 \\leq i,j < n}$ be a stochastic matrix. We define a random process $X_0$, $X_1$, ... as follows. Each $X_t$ is an element of $\\{0, 1, \\ldots, n-1\\}$ (ie the indices of the stochastic matrix are \"states\" of the process). Assuming that $X_t$ is in the state $i$ the variable $X_{t+1}$ is chosen randomly in $\\{0, 1, \\ldots, n-1\\}$ according to the probability vector $(p_{ij})_{0 \\leq j < n}$ (ie the $i$-th row of the matrix $P$).\n", "\n", "Using the function `multi_toss(p)` that you already wrote, write a function `transition(P, i)` where $P$ is a $n \\times n$ stocahstic matrix and $i$ an index in $\\{0, 1, \\ldots, n-1\\}$, returns an integer in $\\{0, 1, \\ldots, n-1\\}$ according to the probabilities given by the $i$-th row of $P$." ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "Check that your function works using the matrices `P1`, `P2`, `P3` and `P4` that were previously defined." ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "Write a function `trajectory(P, x0, t)` that given a $n \\times n$ stochastic matrix, an initial condition $x_0$ and a positive integer $t$ return a simulation of the Markov chain of length $t$ starting from $x_0$." ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "Use your function to get trajectories of length 20 for the four matrices $P_1$, $P_2$, $P_3$ and $P_4$." ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "Make a graphics out of these trajectories (you can use the primitive `point2d` or `line2d`)" ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "## Statistical properties of trajectories\n", "\n", "Given a Markov chain $X_0, X_1, \\ldots$ we will be interested in the number of times the chain visits a given states between time $0$ and $t$. Let us define\n", "\n", "$$S_N(i) := \\# \\{t \\in \\{0, 1, \\ldots, N-1\\}:\\ X_t = i\\}$$\n", "\n", "(Note that $S_n(i)$ is a random variable that depends on the realisation of the Markov chain).\n", "\n", "Write a function `statistics(P, x0, i, N)` that simulates the Markov chain given by $P$ and returns the list of values $S_0(i)$, $S_1(i)$, $S_2(i)$, ..., $S_{N-1}(i)$ obtained along this trajectory." ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "On the same graphics, plot several realisation of the sequence $(S_n(0))_n$ for the Markov chain given by the matrix $P_1$" ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "Repeat the operation with $S_n(0)$ replaced by $S_n(i)$ for all possible values of $i$." ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "What do you observe?\n", "\n", "Make the same experiment with the other matrices $P_2$, $P_3$ and $P_4$." ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "## Equilibrium\n", "\n", "We are now given a stochastic matrix $P$ a fixed initial condition $X_0$ and an integer $t$ and are interested in the random variable $X_t$. That is to say the distribution of the Markov chain after $t$ steps.\n", "\n", "Write a function `distribution(P, x0, t)` that given a matrix $P$, an initial condition $x_0$ and a time $t$ return $x_t$" ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "Write a function `sample(P, x0, t, num)` that given a matrix $P$ an initial condition $x_0$, a time $t$ and a positive integer $num$ return a sample of $num$ simulation of $x_t$" ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "Using samples of size 1000, compute the empirical distribution of $x_100$ for the matrix $P_1$" ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "For $n=0,1,2,\\ldots,100$ compute $u (P_1)^n$ where $u$ is the row vector of floting point numbers $u = (1.0, 0.0, 0.0, 0.0)$" ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "What do you observe?\n", "\n", "Repeat the operation starting from the row vector $u = (0.25, 0.25, 0.25, 0.25)$" ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} }, { "source": [ "Make the same experiment with the matrices $P_2$, $P_3$ and $P_4$" ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# edit here" ], "outputs": [], "metadata": {} } ], "metadata": { "kernelspec": { "display_name": "sagemath", "name": "sagemath" } } }