2016-04-07 1 views
1

Je le code suivant pour une chaîne de Markov:Simuler un chemin de chaîne de Markov 1000 fois

simulation.mc=function(i0,P,n.sim){ 
S=1:nrow(P) 
X=rep(0,n.sim) 
X[1]=i0 
for (n in 2:n.sim){ 
    X[n]=sample(x=S,size=1,prob=P[X[n-1],]) 
    } 
return(X) 
} 

P=matrix(
c(
    0,1/2,0,1/2,0,0,0, 
    1/2,0,1/2,0,0,0,0, 
    0,1,0,0,0,0,0, 
    1/3,0,0,0,1/3,1/3,0, 
    0,0,0,1,0,0,0, 
    0,0,0,1/2,0,0,1/2, 
    0,0,0,0,0,0,1 
),nrow=7,byrow=T);P 

X=simulation.mc(1,P,100) 


T=min(which(X==7)) 

Je dois calculer le nombre moyen d'étapes avant d'atteindre l'état 7.

Je sais que je dois pour exécuter au moins 1000 échantillons du chemin, compter le nombre d'étapes dans chaque échantillon et ensuite calculer la valeur moyenne (bien que certains chemins n'atteignent pas l'état 7).

Je l'ai fait, mais ne fonctionne toujours pas:

n.sim=100 
X[i]=rep(0,n.sim) 
for (i in 1:100) 
{ X[i]=simulation.mc(1,P,100) 
} 

pourquoi cela ne fonctionne pas? Comment puis-je inclure une boucle dans une boucle pour inclure la fonction qui compte le nombre d'étapes os? Merci d'avance pour tout conseil.

Répondre

1

Vous pouvez utiliser replicate au lieu d'une boucle:

replicate(1000, min(which(simulation.mc(1,P,100)==7))) 

@JDB fourni une option pour utiliser une boucle. Voici quelques plus:

# To save each entire chain in a list 
n.sim=100 
X = list() 
for (i in 1:1000) { 
    X[[i]] = simulation.mc(1,P,n.sim) 
} 

# To save just the number of steps to get to 7 
n.sim=100 
X = rep(NA, 1000) 
for (i in 1:1000) { 
    X[i] = min(which(simulation.mc(1,P,n.sim)==7)) 
} 
+0

J'aime les additions la chose intéressante au sujet de la 'min (ce qui())' méthode que je tout à fait comparable est qu'il lance un avertissement si cette simulation ne jamais atteint un 7, ce qui arrive. certains (j'ai simulé combien de fois il se produit dans 1000 répétitions de l'appel original 'replicate()', et a obtenu une moyenne de 13,9 pour 1000.) – JDB

+0

Oui, 'min' retournera' Inf' si la condition n'est jamais remplie. Vous devrez augmenter 'n.sim' pour que le cha dans (presque) toujours atteint 7. J'ai essayé avec 'n.sim = 300' et chaque chaîne a atteint 7 dans 10.000 simulations. – eipi10

+0

Merci beaucoup @ eipi10 et @JDB! Maintenant, je comprends mieux l'utilisation des boucles et 'replicate()'. C'est ainsi que j'ai résolu le problème pour éviter 'Inf' lors du calcul des étapes moyennes:' y <-replicate (1000, min (which (simulation.mc (1, P, 100) == 7))); y2 <-y [est.finite (y)]; moyenne (y2) ' – Maruska

0

Il semble que vous essayez de créer 100 x 100 trame de données, mais vous appelez un vecteur déjà créé à la place.

C'est pourquoi votre appel de X[i]=rep(0,n.sim) ne fonctionne pas. (Je pense que vous pouvez avoir signifié X[1]=rep(0,n.sim), puisque vous définissez seulement i dans la boucle.Vous ne pouvez pas remplir une colonne d'un vecteur de la même manière que vous pouvez un data.frame ...

Essayez :

X <- data.frame(matrix(nrow=100, ncol=100)) 
for (i in 1:100) 
    { X[i]=simulation.mc(1,P,100) 
}