2017-11-16 8 views
1

저는 독자적으로 R을 배우고 있으며, markovchain 패키지를 사용하여 Rstudio에서 전이 확률 행렬을 만드는 데 어려움을 겪고 있습니다. 먼저 DNA 서열의 전이 확률을 계산하려고했습니다.마르코프의 체인 전이 확률 행렬을 만드는 법

ATTCAACACATCCAGCCACATGCTCCGAGAGGAGGCAGAGGGCCCCCGGAATGATGCTTACCGAGATTCTTGTTTTTATCCTCGTGGTTGTTTAAAAACGAGTTGAAACTGACGGCATGTCGGACTATAAGCTACTTACTCACCATAGACGTGACCATAGGCCCTAAAACGTTACCGAGATATTCACTTCTAATAACAGTTGTCGGCAGAGCCAAAAGGCCGGGTGATAATACTTTAAAAAGGGAGTTGATTGTTGTATCTAATCCTAGAATGTCAAGAGCGACCATAACAAGATAATTCGGCAGAGCCAGAAAGCGTTCAAGGACTAGAACCATACCGAGACGCAAACGTTCAGGTCGAACTCTAATACCGATTAGT 

하지만 어떻게 전이 확률 행렬이 같은 순서로 계산 될 수있다, 나는 R 인덱스를 사용하여 생각했지만, 난 정말 그 전이 확률을 계산하는 방법을 모르겠어요.

R에서 이렇게하는 방법이 있습니까? 나는 매트릭스에서 그 확률의 출력이 이런 식으로 뭔가해야한다고 추측하고있다 : DNA

 A T C G 
    A 0.60 0.10 0.10 0.20 
    T 0.10 0.50 0.30 0.10 
    C 0.05 0.20 0.70 0.05 
    G 0.40 0.05 0.05 0.50 

답변

0

만들기

DNA <- "ATTCAACACATCCAGCCACATGCTCCGAGAGGAGGCAGAGGGCCCCCGGAATGATGCTTACCGAGATTCTTGTTTTTATCCTCGTGGTTGTTTAAAAACGAGTTGAAACTGACGGCATGTCGGACTATAAGCTACTTACTCACCATAGACGTGACCATAGGCCCTAAAACGTTACCGAGATATTCACTTCTAATAACAGTTGTCGGCAGAGCCAAAAGGCCGGGTGATAATACTTTAAAAAGGGAGTTGATTGTTGTATCTAATCCTAGAATGTCAAGAGCGACCATAACAAGATAATTCGGCAGAGCCAGAAAGCGTTCAAGGACTAGAACCATACCGAGACGCAAACGTTCAGGTCGAACTCTAATACCGATTAGT" 

문자로 분할 그것은 문자

DNA_list <- unlist(strsplit(DNA, split = "")) 

독특한 요소를 검색

DNA_unique <- unique(DNA_list) 
,

matrix <- matrix(0, ncol = length(DNA_unique), nrow=length(DNA_unique)) 

그것을 기입 행렬 생성 : I ELT를하는 소자 I + 1 매트릭스의 대응하는 셀을 추가.

for (i in 1:(length(DNA_list) - 1)){ 
    index_of_i <- DNA_unique == DNA_list[i] 
    index_of_i_plus_1 <- DNA_unique == DNA_list[i + 1] 
    matrix[index_of_i, index_of_i_plus_1] = matrix[index_of_i, index_of_i_plus_1] + 1 
} 

matrix <- matrix/rowSums(matrix) 

> matrix 
      [,1]  [,2]  [,3]  [,4] 
[1,] 0.3000000 0.2166667 0.2250000 0.2583333 
[2,] 0.3068182 0.2954545 0.2159091 0.1818182 
[3,] 0.2857143 0.2142857 0.2619048 0.2380952 
[4,] 0.3764706 0.2235294 0.1882353 0.2117647 

NB를 정상화 : 당신이 계산하기 위해 정말 많은 DNA가있는 경우 빠른 방법으로 그것을 수행 할 수있는 방법이있을 수 있습니다. 그러나 여기서 그것은 충분히 빠르다고 봅니다.

1

markovchain 패키지를 사용하면 도움이됩니다. 첫째, 데이터가

seq <- "ATTCAACACATCCAGCCACATGCTCCGAGAGGAGGCAGAGGGCCCCCGGAATGATGCTTACCGAGATTCTTGTTTTTATCCTCGTGGTTGTTTAAAAACGAGTTGAAACTGACGGCATGTCGGACTATAAGCTACTTACTCACCATAGACGTGACCATAGGCCCTAAAACGTTACCGAGATATTCACTTCTAATAACAGTTGTCGGCAGAGCCAAAAGGCCGGGTGATAATACTTTAAAAAGGGAGTTGATTGTTGTATCTAATCCTAGAATGTCAAGAGCGACCATAACAAGATAATTCGGCAGAGCCAGAAAGCGTTCAAGGACTAGAACCATACCGAGACGCAAACGTTCAGGTCGAACTCTAATACCGATTAGT" 

그런 다음 패키지

를 사용
library(markovchain) 
base_sequence <- strsplit(seq, "")[[1]] 
mcX <- markovchainFit(base_sequence)$estimate 
mcX 

#   A   C   G   T 
# A 0.3000000 0.2250000 0.2583333 0.2166667 
# C 0.2857143 0.2619048 0.2380952 0.2142857 
# G 0.3764706 0.1882353 0.2117647 0.2235294 
# T 0.3068182 0.2159091 0.1818182 0.2954545 
+0

코드를 언급 할 수 약간 pls는? –