为马尔可夫模型创建一个转移矩阵

Question

为马尔可夫模型创建一个转移矩阵

4

我需要关于马尔可夫链和数据预处理相关主题的帮助。假设我有以下矩阵，它将个体与随时间变化的状态联系起来：

     ID Time1 Time2
1 14021     A     A
2 15031     B     A
3 16452     A     C

我想要获取这个矩阵的状态转移矩阵：因此，所需的是：

同样的事情，但现在根据从该状态转换的总次数进行加权，即

  A    B   C
A 0.5  0  0.5
B 1    0   0
C 0    0   0

（因为有两个从状态A出发的转移）。我知道markovchain包有一个功能，可以在有序列的情况下进行操作，比如AAABBAAABBCC，但是如果数据设置成我这样的形式，就不能使用该功能。理想情况下，直接的程序会很好，但如果有一种方法可以将数据转换为一组序列，那也可以。

- Arrebimbomalho

如果您可以从数据框中创建 c("B","A","A","C")，那么您就可以使用 markovchain::createSequenceMatrix。 - M--

谢谢。但我不太确定如何进行这个序列... - Arrebimbomalho

1

我也是一样。这就是为什么我走了一条长路的原因。 - M--

3个回答

3

使用 igraph 方法，所以使用 Joseph 回答中的 df：

library(igraph)

g <- graph_from_data_frame(df)

E(g)$weight = 1/degree(g, mode="out")[df$Time1] # get counts

as_adj(g, attr = "weight", sparse=FALSE) # output weighted adjacency matrix

    A B   C
A 0.5 0 0.5
B 1.0 0 0.0
C 0.0 0 0.0

- user20650

1

@Arrebimbomalho 这绝对是更好的答案。如果我是你，我会接受这个答案。 - M--

1

@Arrebimbomalho，完全同意Masoud的观点！ - Joseph Wood

1

@JosephWood 我希望你能理解这不是个人攻击。只是认为这个答案更优雅（当然，首先是与我的答案相比）。干杯。 - M--

2

@Masoud，一点也不...我真的是认为这是一个更好的答案。我非常欣赏把社区放在第一位的态度。 - Joseph Wood

2

肯定有更好的方法。这是我在一个无聊的周五下午用循环做的涂鸦。

lvls <- sort(unique(unlist(df[,-1])))

dat <- matrix(0, nrow= length(lvls), ncol= length(lvls))

colnames(dat) <- lvls
rownames(dat) <- lvls

concat <- paste0(df[,2], df[,3])

for (i in 1:length(lvls)) {
  for (j in 1:length(lvls)) {
    dat[i,j] <- paste0(rownames(dat)[i], colnames(dat)[j])
  }
}

dat <- matrix(sapply(dat, function(x) length(grep(x, concat))), 
       nrow= length(lvls), ncol= length(lvls))

colnames(dat) <- lvls
rownames(dat) <- lvls

dat

##   A B C
## A 1 0 1
## B 1 0 0
## C 0 0 0

dat <- dat / rowSums(dat)
dat[is.na(dat)] <- 0

dat

##    A B   C
##A 0.5 0 0.5
##B 1.0 0 0.0
##C 0.0 0 0.0

- M--

太棒了！而且它是可扩展的，所以可能能够处理这个任务（我需要的真实矩阵有数百万条观测数据...） - Arrebimbomalho

@Arrebimbomalho，我在中间有两个循环，可能不会很快，但假设是（1M行*100列），完成起来也不会超过几分钟吧。 - M--

似乎可以使用 table 来完成第一部分，即 df[c("Time1", "Time2")] <- lapply(df[c("Time1", "Time2")], factor, levels=c("A", "B", "C")) ; table(df)。 - user20650

@user20650 尝试创建一个具有 >= 2^31 元素的表格。这只会在处理大型数据集时发生，对吗？ - M--

1

我们只会在这两列上使用“table” - 输出大小将取决于唯一级别的数量。 - user20650

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Joseph Wood · Accepted Answer

这里还有一个与 base R 相关的解决方案。

df <- data.frame(Time1 = c("A","B","A"), Time2 = c("A","A","C"), stringsAsFactors = FALSE)

myStates <- sort(unique(c(df$Time1, df$Time2)))
lenSt <- length(myStates)

currState <- match(df$Time1, myStates)
nextState <- match(df$Time2, myStates)
transMat <- matrix(0L, lenSt, lenSt)

transMat[cbind(currState, nextState)] <- 1L
transMat <- transMat/rowSums(transMat)
transMat[is.na(transMat)] <- 0

transMat
     [,1] [,2] [,3]
[1,]  0.5    0  0.5
[2,]  1.0    0  0.0
[3,]  0.0    0  0.0