P. Carmona
Pour avancer dans le notebook et exécuter les cellules il faut taper Shift+Enter ou utiliser la barre d'outils ci-dessus et choisir Cell, Run Cell and select Below
Vous répondrez aux questions en modifiant ce notebook. En insérant des cellules de type Markdown pour le texte et des cellules de type code pour le code.
Ensuite vous sauvez ce notebook sous le nom Prenom_Nom_tpmarkov.ipynb et vous le déposez sur Moodle
Nous allons simuler un processus de Poisson d'intensité $\lambda=1/10 $ mn sur un intervalle de temps d'une journée $T=24*60=1440$. La commande
import numpy as np
import matplotlib.pyplot as plt
import numpy.random as rnd
%matplotlib inline
lam=0.1
x=rnd.exponential(scale=1/lam)
x
génère une variable exponentielle de paramètre $\lambda$ et la place dans $x$. Pour générer un échantillon de $n=500$ variables, calculer la moyenne et la variance empiriques, et tracer un histogramme normalisé,avec 50 barres, on écrit
n=500
x=rnd.exponential(scale=1/lam,size=n)
x.mean()
x.var()
plt.hist(x, 50, normed=1)
Lorsque l'on désire des informations détaillées sur une fonction on tape
?plt.hist
On améliore le graphique ci-dessus en comparant l'histogramme avec la fonction densité.
plt.hist(x, 50, normed=1)
tt=np.linspace(0.0, 50,200)
plt.plot(tt,lam*np.exp(-lam*tt))
En utilisant la structure de contrôle
while True:
...
if cond:
break
coder une fonction
...
def spoi(lam,T):
qui prend en paramètre $\lambda$ et $T$, et qui renvoie un vecteur des temps de saut du processus de Poisson sur l'intervalle $[0,T]$.
On remarque que
plt.plot(spoi(0.1,100))
n'affiche ansolument pas ce que l'on veut mais la courbe inverse. Réfléchissez bien et affichez la bonne courbe.
Une solution est d'écrire
s=spoi(0.1,100)
plt.step(s,range(len(s)),where='mid')
Coder une seconde fonction sspoi qui a les mêmes paramètres, mais qui utilise la construction générale d'un processus de Poisson. Donc elle genère une variable $N$ de loi $Pois(\lambda T)$ Puis elle génère $N$ uniformes sur $[0,T]$ qu'elle trie. Utiliser les commandes ,rnd.poisson,rnd.runif et x.sort() pour trier un tableau de valeurs. Et affichez également le processus obtenu.
Pour être pleinement convaincu que cette seconde méthode génère bien un processus de Poisson classique, vous allez illustrer le fait que le premier temps de saut est une variable exponentielle de paramètre $\lambda$ en générant $n=500$ telles variables et en traçant l'histogramme et la densité de la loi exponentielle.
Dans un premier temps on suppose que l'on arrive à $t=15h00=15*60$ minutes et que l'on désire déterminer les valeurs moyennes de l'âge $A_t=t-S_{N_t}$ et du temps de vie résiduel $R_t = S_{1+N(t)}-t$, qui sont respectivement l'intervalle depuis le passage du dernier bus et le temps d'attente du prochain bus.
On écrit donc une fonction attente qui prend en paramètres $\lambda$ et $t$ et revoie le couple de valeurs $A_t,R_t$.
def attente(lam,t):
...
return ((at,rt))
Puis on utilise
lam=0.1
t=300
res=[attente(lam,t) for i in range(1000)]
qui fabrique un vecteur contenant 1000 réplications (i.e. les valeurs de 1000 appels successifs de la fonction) et on les place dans une matrice, en calculant la moyenne des lignes ce qui par la loi des grands nombres donne des estimations des moyennes recherchées
m=res.reshape(2,500)
m.mean(axis=1)
Pour illustrer l'indépendance par rapport à l'instant d'arrivée, faire les mêmes calculs en tirant cette fois au hasard votre temps d'arrivée à l'arrêt, en utilisant par exemple une loi uniforme sur 15h, 15h15.
Vérifier vos simulations en utilisant le calcul matriciel. Sous Python, le produit matriciel des matrices p et q s'écrit np.dot(p,q) ce qui est très différent de $p*q$. Ecrivez une fonction puis(p,n) qui renvoie la matrice $p^n$. (Vous pouvez utiliser un algorithme de calcul rapide de puissance si vous savez le coder).
On fixe les valeurs du coût de vidage $K=100$ euros et $\lambda=1$, i.e. on a un message en moyenne par unité de temps. En faisant varier le coût unitaire $h$ de stockage d'un message (valeurs $0.1,1,10 $ ou $100$), calculer le coût moyen par cycle en effectuant un grand nombre (n=500) de cycles, et ce pour suffisamment de valeurs de $T$ pour identifier la valeur optimale.
On pourra écrire d'abord une fonction
def coutcycle(T,h,K,lam):
qui calcule le coût unitaire observé pour la réalisation d'un cycle.
Commandes Python : np.linspace pour construire un vecteurs de valeurs uniformément réparties
Ensuite il faut une fonction qui donne le coût moyen estime par la loi forte des grands nombres
def coutmoyenempi(T,h,K,lam,n):
Puis on tracera, au moins pour une valeur de $h$, la courbe qui donne le coût moyen estimé en fonction de $T$. Et la on vérifiera que la valeur $T^*$ calculée en TD est bien optimale.