如何解决是否有 R 包可以从频率表中计算一阶转移矩阵?
我有一个由 8 亿条记录聚合而成的频率表,我想知道是否可以使用包从频率表中计算一阶转移矩阵,这是不对称的,因为某些状态再也没有发生过。频率表的一个样本是:
library(data.table)
model.data <- data.table(state1 = c(3,1,2,3),state2 = c(1,2),Freq = c(1,3,4))
model.data 看起来像这样:
state1 | state2 | n |
---|---|---|
3 | 1 | 1 |
1 | 2 | 2 |
2 | 1 | 3 |
3 | 2 | 4 |
使用包 pollster,我可以计算出比例表:
library(pollster)
crosstab(model.data,state1,state2,Freq)
state1 | 1 | 2 | n |
---|---|---|---|
1 | 0 | 100 | 2 |
2 | 100 | 0 | 3 |
3 | 20 | 80 | 5 |
然而,我正在寻找的对称转移矩阵是:
state1 | 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 表示。 (但也许不是因为您发布了增加了额外混淆列的民意调查代码。)在这里,我们将反转该过程,以便我们可以恢复矩阵(R table
对象继承自该矩阵)。
请注意,我使用的是您的数据对象,但没有使用 pkg:pollster 代码,因为您的表似乎不是基于该 data.table 对象。
如何获得零列,...只需在 state2=3
“列”位置放入单个零数据元素。您只需要为整个列在 state2 中添加一个数据点,但显然它需要来自某个 state1 值。它可以来自任何状态 1 值:
model.data <- data.table(state1 = c(3,1,2,3,3),state2 = c(1,Freq = c(1,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
补充说明:只是为了表明这在“民意调查”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 个包。 (他们似乎假设每个人都知道该怎么做,尽管当我为上面链接的答案编写代码时,我需要回到我的书上复习。)
带有 igraph
的选项
model.data %>%
setorder(state1) %>%
graph_from_data_frame() %>%
as_adjacency_matrix(attr = "Freq",sparse = FALSE) %>%
proportions(1) # 1 sets rows as the margin,similar to `prop.table`
给予
1 2 3
1 0.0 1.0 0
2 1.0 0.0 0
3 0.2 0.8 0
或使用基数 R
> proportions(xtabs(Freq ~ .,model.data),1)
state2
state1 1 2
1 0.0 1.0
2 1.0 0.0
3 0.2 0.8
版权声明:本文内容由互联网用户自发贡献,该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 dio@foxmail.com 举报,一经查实,本站将立刻删除。