Simulation d'une CMTD simple

set.seed(1)

On simule le nombre de clients dans le système à l'aide des équations d'évolution vues en cours.

nb_clients = function(T=1000,a=0.5,p=0.95){ #T=horizon temporel / nb d'événements; a=proba d'arrivée; p=proba de transmission réussie

  #etat de la simulation
  X<-c() #trajectoire (nombre de clients)
  A<-c() #instants d'arrivées 
  D<-c() #instants de departs 

  A=rbinom(T,1,a)
  D=rbinom(T,1,p) #en cas de départ d'une file vide, il ne se passera rien
  X[1]<-0 #état initial
  for(n in 2:T){  #n est le temps simulé
    X[n]=X[n-1]+A[n]-min(D[n],X[n-1]) #fonction de transition Phi
  }
  X
}

Exemple de trajectoire. On prend un horizon pas trop long pour visualiser les transitions.

T<-100
a<-0.5
p<-a+0.005
X=nb_clients(T,a,p)
plot(seq(1:T),X,xlab="itérations",ylab="Nombre de clients dans le système",type="s")

plot of chunk unnamed-chunk-3

La file ne se remplit pas beaucoup (malgré le faible écart entre p et a).

Essayons de voir si la file “n'explose pas” au bout d'un certain temps (plus long que 100!) avec une file stable mais chargée (a/p proche de 1).

T=10000
a=0.9
p=a+0.0005
X=nb_clients(T,a,p)
mean(X)
## [1] 14.3029
max(X)
## [1] 38

On voudrait aussi connaître la distribution du nombre de clients (attention, ici sur 1 seule trajectoire !) :

hist(X,main="Répartition des états rencontrés par la trajectoire", xlab="états")

plot of chunk unnamed-chunk-5 Difficile à interpréter…

On voudrait regarder le comportement sur plusieurs trajectoires.

ech=c()
T<-10000
R<-100 #nb de répétitions
for(i in seq(1:R)){
  ech<-c(ech,mean(nb_clients(T,a,p)))
}
hist(ech,main="Répartition du taux d'occupation moyen", xlab="E[X]")

plot of chunk unnamed-chunk-6