2010-05-02 2 views
8

Ищу что-то вроде «MSM» пакета, но для дискретных цепей Маркова. Например, если бы у меня была матрица перехода, определенная как таковаяR библиотека для дискретных цепей Маркова моделирования

Pi <- matrix(c(1/3,1/3,1/3, 
0,2/3,1/6, 
2/3,0,1/2)) 

для состояний A, B, C. Как я могу имитировать цепь Маркова в соответствии с этой матрицей перехода?

Спасибо,

+2

NVM: нашел ответ http://books.google.com/books?id=AALhk_mt7SYC&lpg=PA116&dq=r%20markov%20chain&pg=PA119#v=onepage&q=r%20markov%20chain&f=false – stevejb

ответ

8

Некоторое время назад я написал набор функций для моделирования и оценки дискретных марковских матриц вероятностей цепи: http://www.feferraz.net/files/lista/DTMC.R.

Соответствующий код для того, что вы спрашиваете:

simula <- function(trans,N) { 
     transita <- function(char,trans) { 
       sample(colnames(trans),1,prob=trans[char,]) 
     } 

sim <- character(N) 
sim[1] <- sample(colnames(trans),1) 
for (i in 2:N) { 
    sim[i] <- transita(sim[i-1],trans) 
} 

sim 
} 

#example 
#Obs: works for N >= 2 only. For higher order matrices just define an 
#appropriate mattrans 
mattrans <- matrix(c(0.97,0.03,0.01,0.99),ncol=2,byrow=TRUE) 
colnames(mattrans) <- c('0','1') 
row.names(mattrans) <- c('0','1') 
instancia <- simula(mattrans,255) # simulates 255 steps in the process 
6

Argh, вы нашли решение в то время как я пишу это для вас. Вот простой пример, который я придумал:

run = function() 
{ 
    # The probability transition matrix 
    trans = matrix(c(1/3,1/3,1/3, 
       0,2/3,1/3, 
       2/3,0,1/3), ncol=3, byrow=TRUE); 

    # The state that we're starting in 
    state = ceiling(3 * runif(1, 0, 1)); 
    cat("Starting state:", state, "\n"); 

    # Make twenty steps through the markov chain 
    for (i in 1:20) 
    { 
     p = 0; 
     u = runif(1, 0, 1); 

     cat("> Dist:", paste(round(c(trans[state,]), 2)), "\n"); 
     cat("> Prob:", u, "\n"); 

     newState = state; 
     for (j in 1:ncol(trans)) 
     { 
      p = p + trans[state, j]; 
      if (p >= u) 
      { 
       newState = j; 
       break; 
      } 
     } 

     cat("*", state, "->", newState, "\n"); 
     state = newState; 
    } 
} 

run(); 

Обратите внимание, что ваша матрица вероятности перехода не подводить к 1 в каждой строке, что он должен делать. Мой пример имеет слегка измененную матрицу перехода вероятности, которая придерживается этого правила.

+0

Спасибо за ответ. Ваш код очень читабельен. Я очень ценю это. – stevejb

+2

Читаемый код? По моему опыту, эта концепция полностью потеряна для большинства людей, которые используют «R». Надеюсь, поможет! – icio

+4

Для того, чтобы сгенерировать случайное число от 1 до 3, по-моему 'образец (1: 3, 1)' немного легче, чем 'обращал внимания потолку (3 * runif (1, 0, 1))'. Кроме того, для самого внутреннего цикла for вы можете просто использовать 'newState <- sample (1: ncol (trans), 1, prob = trans [state,])'. Это более ясно показывает, что происходит. И тогда вам даже не придется нормализовать строки. –

5

Вы можете теперь использовать markovchain пакет, доступный в CRAN. Пользователь manual. довольно неплохо и имеет несколько примеров.

Смежные вопросы