Commit 6b974123 authored by vincentvigon's avatar vincentvigon
Browse files

3 premier chap de 'alea'

parent a005fb67
Pipeline #2564 failed with stage
in 2 minutes and 59 seconds
......@@ -3,4 +3,4 @@
**/__pycache__
**/.ipynb_checkpoints/
**/.DS_Store
**/build
\ No newline at end of file
......@@ -5,8 +5,10 @@
<inspection_tool class="PyCompatibilityInspection" enabled="true" level="WARNING" enabled_by_default="true">
<option name="ourVersions">
<value>
<list size="1">
<list size="3">
<item index="0" class="java.lang.String" itemvalue="3.5" />
<item index="1" class="java.lang.String" itemvalue="3.6" />
<item index="2" class="java.lang.String" itemvalue="3.7" />
</list>
</value>
</option>
......
......@@ -4,4 +4,7 @@
<option name="languageLevel" value="ES6" />
</component>
<component name="ProjectRootManager" version="2" project-jdk-name="Python 3.6" project-jdk-type="Python SDK" />
<component name="PythonCompatibilityInspectionAdvertiser">
<option name="version" value="3" />
</component>
</project>
\ No newline at end of file
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
......@@ -19,302 +19,6 @@
"np.set_printoptions(precision=3,suppress=True)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Préliminaires d'algèbre linéaire\n",
"\n",
"Attention, en numpy les vecteurs, les matrices lignes et les matrices colonnes sont des objets différents:\n",
"\n",
"* vecteur.shape = (?)\n",
"* matrice_ligne.shape = (1,?)\n",
"* matrice_colonne.shape = (?,1)\n",
"* matrice_quelconque.shape = (?,?)\n",
"\n",
"Observez bien les sorties consoles: les vecteurs s'écrivent avec 1 crochet, les matrices avec 2 crochets.\n",
"\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Multiplication matricielle\n",
"\n",
"`np.matmul()` s'applique uniquement entre matrices."
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"mat_col\n",
" [[1.]\n",
" [1.]\n",
" [1.]]\n",
"\n",
"mat_lin\n",
" [[1. 1. 1.]]\n",
"\n",
"mat_square\n",
" [[1. 1. 1.]\n",
" [1. 1. 1.]\n",
" [1. 1. 1.]]\n",
"\n",
"mat_square × mat_col\n",
" [[3.]\n",
" [3.]\n",
" [3.]]\n",
"\n",
"mat_lin × mat_square\n",
" [[3. 3. 3.]]\n",
"\n",
"mat_lin × mat_square × mat_col\n",
" [[9.]]\n"
]
}
],
"source": [
"size=3\n",
"mat_col=np.ones(shape=[size,1])\n",
"mat_lin=np.ones(shape=[1,size])\n",
"mat_square=np.ones(shape=[size,size])\n",
"\n",
"print(\"\\nmat_col\\n\",mat_col)\n",
"print(\"\\nmat_lin\\n\",mat_lin)\n",
"print(\"\\nmat_square\\n\",mat_square)\n",
"\n",
"print(\"\\nmat_square × mat_col\\n\",np.matmul(mat_square,mat_col) )\n",
"print(\"\\nmat_lin × mat_square\\n\",np.matmul(mat_lin,mat_square))\n",
"print(\"\\nmat_lin × mat_square × mat_col\\n\", np.matmul(np.matmul(mat_lin,mat_square),mat_col))\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"sinon on peut utiliser `np.dot()` qui permet les multiplications matrice × vecteur, matrice × matrice, vecteur × vecteur (=produit scalaire)"
]
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"vec0 × vec1\n",
" 6.0\n",
"\n",
"vec0 × mat_square\n",
" [3. 3. 3.]\n",
"\n",
"mat_square × vec0\n",
" [3. 3. 3.]\n",
"\n",
"mat_square × mat_square\n",
" [[3. 3. 3.]\n",
" [3. 3. 3.]\n",
" [3. 3. 3.]]\n"
]
}
],
"source": [
"size=3\n",
"vec0=np.ones(shape=[size])\n",
"vec1=2*np.ones(shape=[size])\n",
"mat_square=np.ones(shape=[size,size])\n",
"\n",
"print(\"\\nvec0 × vec1\\n\",np.dot(vec0,vec1))\n",
"print(\"\\nvec0 × mat_square\\n\",np.dot(vec0,mat_square))\n",
"print(\"\\nmat_square × vec0\\n\",np.dot(mat_square,vec0))\n",
"print(\"\\nmat_square × mat_square\\n\",np.dot(mat_square,mat_square))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"et si vous avez python 3.5+, vous pouver utilisez l'opérateur @ qui rend les codes plus lisibles"
]
},
{
"cell_type": "code",
"execution_count": 18,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"vec0 × vec1\n",
" 6.0\n",
"\n",
"vec0 × mat_square\n",
" [3. 3. 3.]\n",
"\n",
"mat_square × vec0\n",
" [3. 3. 3.]\n",
"\n",
"mat_square × mat_square\n",
" [[3. 3. 3.]\n",
" [3. 3. 3.]\n",
" [3. 3. 3.]]\n"
]
}
],
"source": [
"print(\"\\nvec0 × vec1\\n\",vec0 @ vec1)\n",
"print(\"\\nvec0 × mat_square\\n\",vec0 @ mat_square)\n",
"print(\"\\nmat_square × vec0\\n\",mat_square @ vec0)\n",
"print(\"\\nmat_square × mat_square\\n\",mat_square @ mat_square)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"***Exo:*** Multipliez matriciellement des matrices de taille non-compatible.\n",
"Extrayez la partie intéressante du message d'erreur."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### inverse et transposée \n",
"\n",
"Notez que le pseudo-inverse permet d'inverser les matrices non-inversibles (testez). "
]
},
{
"cell_type": "code",
"execution_count": 19,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"mat\n",
" [[1. 0. 2.]\n",
" [0. 1. 0.]\n",
" [0. 0. 1.]]\n",
"\n",
"mat^T\n",
" [[1. 0. 0.]\n",
" [0. 1. 0.]\n",
" [2. 0. 1.]]\n",
"\n",
"mat^(-1)\n",
" [[ 1. 0. -2.]\n",
" [ 0. 1. 0.]\n",
" [ 0. 0. 1.]]\n",
"\n",
"mat^(-1)-pseudo-inverse\n",
" [[ 1. 0. -2.]\n",
" [ 0. 1. 0.]\n",
" [ 0. 0. 1.]]\n"
]
}
],
"source": [
"size=3\n",
"mat=np.zeros(shape=[size,size])\n",
"for i in range(size):\n",
" mat[i,i]=1\n",
"mat[0,size-1]=2\n",
"\n",
"\n",
"print(\"\\nmat\\n\",mat )\n",
"print(\"\\nmat^T\\n\",mat.T)\n",
"print(\"\\nmat^(-1)\\n\",np.linalg.inv(mat))\n",
"print(\"\\nmat^(-1)-pseudo-inverse\\n\",np.linalg.pinv(mat))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### vecteur propre et valeur propre\n",
"\n",
"Souvenez-vous qu'en anglais valeur/vecteur 'propre' c'est 'eigen' value/vector"
]
},
{
"cell_type": "code",
"execution_count": 20,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"mat\n",
" [[0.085 0.287 0.206 0.422]\n",
" [0.191 0.262 0.253 0.295]\n",
" [0.244 0.07 0.06 0.626]\n",
" [0.383 0.001 0.161 0.454]]\n",
"\n",
"eigVal\n",
" [ 1. 0.122 -0.195 -0.066]\n",
"\n",
"eigVec\n",
" [[-0.5 -0.416 -0.418 0.39 ]\n",
" [-0.5 -0.763 0.307 0.432]\n",
" [-0.5 0.402 -0.739 -0.813]\n",
" [-0.5 0.287 0.43 -0.036]]\n"
]
}
],
"source": [
"size=4\n",
"\"\"\"une matrice aléatoire auquel on fait subir une opération de normalisation (laquelle?)\"\"\"\n",
"mat=np.random.uniform(low=0,high=10,size=[size,size])\n",
"rowSum=np.sum(mat,axis=1)\n",
"\"\"\" mat[i,j]=mat[i,j]/rowSum(i) \"\"\"\n",
"for i in range(size): mat[i,:]/=rowSum[i]\n",
"\n",
"eigVal,eigVec=np.linalg.eig(mat)\n",
"print(\"\\nmat\\n\", mat )\n",
"print(\"\\neigVal\\n\",np.real(eigVal))\n",
"print(\"\\neigVec\\n\", np.real(eigVec))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"***A vous:*** \n",
"\n",
"* Pourquoi 1 est-il valeur propre ? Quel est le vecteur propre associé? `np.linalg.eig` donne-t-il les vecteurs propres à droite ou à gauche? Sont-ils disposés en ligne ou bien en colonne?\n",
"\n",
"* Toutes les matrices sont-elles diagonalisables? A quoi servent les np.real()? \n",
"\n",
"* laquelle de ces deux formules est la bonne ?\n",
"\n",
" Diag = P^-1 @ mat @ P \n",
"\n",
"ou\n",
"\n",
" Diag = P @ mat @ P^-1\n",
" \n",
"Vérifiez avec python. Aide: Pour transformez un vecteur en matrice diagonale, utilisez `np.diag(vecteur)`"
]
},
{
"cell_type": "markdown",
"metadata": {},
......
%% Cell type:markdown id: tags:
# Chaine de Markov
%% Cell type:code id: tags:
``` python
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
np.set_printoptions(precision=3,suppress=True)
```
%% Cell type:markdown id: tags:
## Préliminaires d'algèbre linéaire
Attention, en numpy les vecteurs, les matrices lignes et les matrices colonnes sont des objets différents:
* vecteur.shape = (?)
* matrice_ligne.shape = (1,?)
* matrice_colonne.shape = (?,1)
* matrice_quelconque.shape = (?,?)
Observez bien les sorties consoles: les vecteurs s'écrivent avec 1 crochet, les matrices avec 2 crochets.
%% Cell type:markdown id: tags:
### Multiplication matricielle
`np.matmul()` s'applique uniquement entre matrices.
%% Cell type:code id: tags:
``` python
size=3
mat_col=np.ones(shape=[size,1])
mat_lin=np.ones(shape=[1,size])
mat_square=np.ones(shape=[size,size])
print("\nmat_col\n",mat_col)
print("\nmat_lin\n",mat_lin)
print("\nmat_square\n",mat_square)
print("\nmat_square × mat_col\n",np.matmul(mat_square,mat_col) )
print("\nmat_lin × mat_square\n",np.matmul(mat_lin,mat_square))
print("\nmat_lin × mat_square × mat_col\n", np.matmul(np.matmul(mat_lin,mat_square),mat_col))
```
%%%% Output: stream
mat_col
[[1.]
[1.]
[1.]]
mat_lin
[[1. 1. 1.]]
mat_square
[[1. 1. 1.]
[1. 1. 1.]
[1. 1. 1.]]
mat_square × mat_col
[[3.]
[3.]
[3.]]
mat_lin × mat_square
[[3. 3. 3.]]
mat_lin × mat_square × mat_col
[[9.]]
%% Cell type:markdown id: tags:
sinon on peut utiliser `np.dot()` qui permet les multiplications matrice × vecteur, matrice × matrice, vecteur × vecteur (=produit scalaire)
%% Cell type:code id: tags:
``` python
size=3
vec0=np.ones(shape=[size])
vec1=2*np.ones(shape=[size])
mat_square=np.ones(shape=[size,size])
print("\nvec0 × vec1\n",np.dot(vec0,vec1))
print("\nvec0 × mat_square\n",np.dot(vec0,mat_square))
print("\nmat_square × vec0\n",np.dot(mat_square,vec0))
print("\nmat_square × mat_square\n",np.dot(mat_square,mat_square))
```
%%%% Output: stream
vec0 × vec1
6.0
vec0 × mat_square
[3. 3. 3.]
mat_square × vec0
[3. 3. 3.]
mat_square × mat_square
[[3. 3. 3.]
[3. 3. 3.]
[3. 3. 3.]]
%% Cell type:markdown id: tags:
et si vous avez python 3.5+, vous pouver utilisez l'opérateur @ qui rend les codes plus lisibles
%% Cell type:code id: tags:
``` python
print("\nvec0 × vec1\n",vec0 @ vec1)
print("\nvec0 × mat_square\n",vec0 @ mat_square)
print("\nmat_square × vec0\n",mat_square @ vec0)
print("\nmat_square × mat_square\n",mat_square @ mat_square)
```
%%%% Output: stream
vec0 × vec1
6.0
vec0 × mat_square
[3. 3. 3.]
mat_square × vec0
[3. 3. 3.]
mat_square × mat_square
[[3. 3. 3.]
[3. 3. 3.]
[3. 3. 3.]]
%% Cell type:markdown id: tags:
***Exo:*** Multipliez matriciellement des matrices de taille non-compatible.
Extrayez la partie intéressante du message d'erreur.
%% Cell type:markdown id: tags:
### inverse et transposée
Notez que le pseudo-inverse permet d'inverser les matrices non-inversibles (testez).
%% Cell type:code id: tags:
``` python
size=3
mat=np.zeros(shape=[size,size])
for i in range(size):
mat[i,i]=1
mat[0,size-1]=2
print("\nmat\n",mat )
print("\nmat^T\n",mat.T)
print("\nmat^(-1)\n",np.linalg.inv(mat))
print("\nmat^(-1)-pseudo-inverse\n",np.linalg.pinv(mat))
```
%%%% Output: stream
mat
[[1. 0. 2.]
[0. 1. 0.]
[0. 0. 1.]]
mat^T
[[1. 0. 0.]
[0. 1. 0.]
[2. 0. 1.]]
mat^(-1)
[[ 1. 0. -2.]
[ 0. 1. 0.]
[ 0. 0. 1.]]
mat^(-1)-pseudo-inverse
[[ 1. 0. -2.]
[ 0. 1. 0.]
[ 0. 0. 1.]]
%% Cell type:markdown id: tags:
### vecteur propre et valeur propre
Souvenez-vous qu'en anglais valeur/vecteur 'propre' c'est 'eigen' value/vector
%% Cell type:code id: tags:
``` python
size=4
"""une matrice aléatoire auquel on fait subir une opération de normalisation (laquelle?)"""
mat=np.random.uniform(low=0,high=10,size=[size,size])
rowSum=np.sum(mat,axis=1)
""" mat[i,j]=mat[i,j]/rowSum(i) """
for i in range(size): mat[i,:]/=rowSum[i]
eigVal,eigVec=np.linalg.eig(mat)
print("\nmat\n", mat )
print("\neigVal\n",np.real(eigVal))
print("\neigVec\n", np.real(eigVec))
```
%%%% Output: stream
mat
[[0.085 0.287 0.206 0.422]
[0.191 0.262 0.253 0.295]
[0.244 0.07 0.06 0.626]
[0.383 0.001 0.161 0.454]]
eigVal
[ 1. 0.122 -0.195 -0.066]
eigVec
[[-0.5 -0.416 -0.418 0.39 ]
[-0.5 -0.763 0.307 0.432]
[-0.5 0.402 -0.739 -0.813]
[-0.5 0.287 0.43 -0.036]]
%% Cell type:markdown id: tags:
***A vous:***
* Pourquoi 1 est-il valeur propre ? Quel est le vecteur propre associé? `np.linalg.eig` donne-t-il les vecteurs propres à droite ou à gauche? Sont-ils disposés en ligne ou bien en colonne?
* Toutes les matrices sont-elles diagonalisables? A quoi servent les np.real()?
* laquelle de ces deux formules est la bonne ?
Diag = P^-1 @ mat @ P
ou
Diag = P @ mat @ P^-1
Vérifiez avec python. Aide: Pour transformez un vecteur en matrice diagonale, utilisez `np.diag(vecteur)`
%% Cell type:markdown id: tags:
Attention, lorsqu'on part d'une matrice symétrique (ou hermitienne), pour obtenir la diagonalisation dans une base othornormale il faut utiliser `np.linalg.eigh` (h-comme hermitienne). Vérifiez-le!
%% Cell type:markdown id: tags:
## Définition intuitive d'une chaine de Markov
%% Cell type:markdown id: tags:
***NB:*** Le mot "chaine de Markov" (tout court) est pour nous synonime de "chaine de Markov homogène en temps".
On considère un ensemble dénombrable $E$ : l'ensemble des "états". Par exemple $E=\{0,1,2,3,4,5\}$. On considère aussi que les éléments de $E$ sont des sommets d'un graphe, dont les flèches sont pondérées.
![graph_ponder](./img/graph_ponder.png)
Une chaine de Markov $t\to X_t$ est un processus aléatoire qui se ballade en suivant le graphe. Au temps $t=0$ elle se trouve en une sommet $X_0$ donné. Si au temps $t$ elle se trouve en $X_t$, elle tire une des flèches partant de $X_t$, avec une probabilité proportionnelle aux pondérations, et elle suit cette flèche pour arriver à un nouvel état $X_{t+1}$
Notations:
* Le poids mis sur la fléche allant de $x$ à $y$ est noté $P(x,y)$.
* Quand il n'y a pas de flèche entre $x$ et $y$, on note $P(x,y)=0$.
* La matrice $\big ( P(x,y) \big)_{x,y \in E}$ est appelée matrice de transition.
***Remarques:*** Ici, les étiquettes que l'on a mises sur les états n'ont aucune importance. On aurait pu tout aussi bien prendre $E=\{A,B,C,E,D,F\}$. Dans les exemples suivants, le fait que les indexes soient des entiers aura son importance.
%% Cell type:code id: tags:
``` python
def premier_chaine():
P=np.zeros([6,6])
P[0,3]=4
P[1,0]=2.1
P[2,1]=2.5
P[3,1]=3
P[3,4]=2
P[4,1]=2
P[4,2]=0.2
P[4,5]=6
P[5,2]=7.3
sumLine=np.sum(P,axis=1)
""" P[i,j]=P[i,j]/sumLine[i] """
P/=np.expand_dims(sumLine,axis=1)
return P
print(premier_chaine())
```
%%%% Output: stream
[[0. 0. 0. 1. 0. 0. ]
[1. 0. 0. 0. 0. 0. ]
[0. 1. 0. 0. 0. 0. ]
[0. 0.6 0. 0. 0.4 0. ]
[0. 0.244 0.024 0. 0. 0.732]
[0. 0. 1. 0. 0. 0. ]]
%% Cell type:code id: tags:
``` python
def markov_from_P(t_max,P,x0):
X=np.zeros(t_max,dtype=np.int)
X[0]=x0
for t in range(t_max-1):
X[t+1]=np.random.choice(a=range(len(P)),p=P[X[t],:])
return X
t_max=300
P=premier_chaine()
X=markov_from_P(t_max,P,3)
plt.figure(figsize=(20,3))
plt.plot(range(t_max),X);
plt.grid()
```
%%%% Output: display_data