我有一个从 8 亿条记录中聚合而成的频率表,我想知道是否可以使用一个包来计算频率表中的一阶转移矩阵,该矩阵不是对称的,因为某些状态再也没有发生过。频率表的示例是:
library(data.table)
model.data <- data.table(state1 = c(3, 1, 2, 3), state2 = c(1, 2, 1, 2), Freq = c(1,2,3,4))
model.data 如下所示:
状态1 | 状态2 | n |
---|---|---|
3 | 1 | 1 |
1 | 2 | 2 |
2 | 1 | 3 |
3 | 2 | 4 |
使用包 pollster,我可以计算比例表:
library(pollster)
crosstab(model.data, state1, state2, Freq)
状态1 | 1 | 2 | n |
---|---|---|---|
1 | 0 | 100 | 2 |
2 | 100 | 0 | 3 |
3 | 20 | 80 | 5 |
但是,我正在寻找的对称转移矩阵是:
状态1 | 1 | 2 | 3 | n |
---|---|---|---|---|
1 | 0 | 100 | 0 | 2 |
2 | 100 | 0 | 0 | 3 |
3 | 20 | 80 | 0 | 5 |
也就是说,即使没有人转换到状态 3,我仍然想包含状态 3,并且代码应该能够自动找出 3 需要附加一列 0。
由于内存限制和计算速度慢,我不确定带有 markovchainFit 函数的 markovchain 包是否能够处理我需要转换为数百万个序列列表的 8 亿行数据。
有人知道吗?
最佳答案
看来您可能已经知道 stats::xtabs
函数,因为您要求我们处理的结果似乎是 base::as 的结果。 data.frame.table
函数将 table
调用的“宽”结果转换为相同数据的“长”data.frame 表示。 (但也许不是,因为您发布的 pollster 代码添加了一个额外的令人困惑的列。)在这里,我们将反转该过程,以便我们可以恢复一个矩阵(R table
对象继承自该矩阵)。
Do notice that I'm using your data object, but not using pkg:pollster code, since your tables didn't appear to be based on that data.table object.
如何获得零列,...只需在 state2=3
“列”位置放入单个零数据元素。您只需要在 state2 中为整列添加一个数据点,但它显然需要来自某个 state1 值。它可以来自任何状态 1 值:
model.data <- data.table(state1 = c(3, 1, 2, 3, 3),
state2 = c(1, 2, 1, 2, 3),
Freq = c(1,2,3,4, 0))
xtabs(Freq~state1+state2, model.data)
#------------
state2
state1 1 2 3
1 0 2 0
2 3 0 0
3 1 4 0
添加注释:只是为了表明这在“pollster”tidyverse 环境中有效......
> library(pollster)
> crosstab(model.data, state1, state2, Freq)
# A tibble: 3 x 5
state1 `1` `2` `3` n
<dbl> <dbl> <dbl> <dbl> <dbl>
1 1 0 100 0 2
2 2 100 0 0 3
3 3 20 80 0 5
进一步注意,如果您想制作一个转换矩阵,则需要删除“n”列。 (我不太清楚它代表什么。)
关于如何制作转换矩阵(如果需要,则将矩阵除以 rowSums
结果,因为转换矩阵需要使每行总和为单位)
mat <- xtabs(Freq~state1+state2, model.data)
trans_mat <- mat/rowSums(mat)
trans_mat
#-----
state2
state1 1 2 3
1 0.0 1.0 0.0
2 1.0 0.0 0.0
3 0.2 0.8 0.0
现在您可以使用矩阵乘法计算任何离散间隔的状态:参见 ?'%*%'
或矩阵求幂 ?expm::expm
这里进一步编码了一个与转移矩阵上的矩阵运算相关的图,以生成马尔可夫模拟: Simple Markov Chain in R (visualization)
markovchain
包中提供了对马尔可夫序列的进一步统计操作,但我没有看到它有任何实际从数据构造转换矩阵的功能。我可能错了,因为我只阅读了小插图的前 5 个包。 (他们似乎认为每个人都知道如何做到这一点,尽管当我为上面链接的答案编写代码时,我需要回到我的书本上进行复习。)
关于r - 是否有一个 R 包可以从频率表中计算一阶转换矩阵?,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/68603750/