R中计算平稳分布的源代码答案

【问题标题】：Source code for calculation of stationary distribution in RR中计算平稳分布的源代码
【发布时间】：2019-05-10 00:06:15
【问题描述】：

看看this link。

我正在尝试理解以下用于查找矩阵的平稳分布的源代码：

# Stationary distribution of discrete-time Markov chain
# (uses eigenvectors)
stationary <- function(mat)
{
    x = eigen(t(mat))$vectors[,1]
    as.double(x/sum(x))
}

我自己测试了以下源代码：

> rm(list=ls())
>
> P <- matrix(c(0.66, 0.34,
+               0.66, 0.34), nrow=2, ncol=2, byrow = TRUE)
>
> x <- eigen(t(P))
> x$values
[1] 1 0

$vectors
          [,1]       [,2]
[1,] 0.8889746 -0.7071068
[2,] 0.4579566  0.7071068

> y <- x$vectors[,1]
> y
[1] 0.8889746 0.4579566
>

看起来像命令

y <- x$vectors[,1]

正在选择矩阵的第一列。

为什么不简单地写成下面这样？

# Stationary distribution of discrete-time Markov chain
# (uses eigenvectors)
stationary <- function(mat)
{
    x = eigen(t(mat))
    y = x[,1]
    as.double(y/sum(y))
}

引入美元符号和向量关键字的原因是什么？

【问题讨论】：

美元运算符继承自 S/Splus。最好避免它。（你的观察是正确的，IMO）
“引入美元符号和矢量关键字的原因是什么？” 我很困惑你在问什么。 eigen 返回一个list，其中元素vectors 是存储在矩阵中的特征向量。因此，您只需从eigen 返回的list 访问vectors 元素。使用$ 索引访问list 元素没有任何问题，例如lst <- list(a=1, b=2); lst$a.
@user366312 我认为 wildplasser 可能误解了您的问题。正如我所说，使用$ 访问来自命名list 的元素绝对没有任何问题。你也会使用$-indexing 来访问来自data.frame 的元素（列），不是吗？
嗯。也许我误解了这个问题。使用$ 表示法访问命名list 元素是R 中的绝对标准做法。有关涉及eigen 的具体示例，请参见例如Eigenvalues and eigenvectors from the official "An Introduction to R" on CRAN; $ 索引当然不在 R 中应该避免。
我认为您提出的更简单的代码不会产生正确的输出。如前所述，eigen 返回一个包含两个对象的列表。你必须告诉 R 你想要列表中的哪个项目（这里是vectors）然后告诉 R 你想要矩阵的哪一部分。 x[[2]][,1] 将等同于 x$vectors[,1] 。这两种都很常见。

标签： r markov-chains markov

【解决方案1】：

让我们测试一下你的提议：

> P <- matrix(c(0.66, 0.34, 0.66, 0.34), nrow=2, ncol=2, byrow = TRUE)
> x <- eigen(t(P))
> print(x)
eigen() decomposition
$values
[1] 1 0

$vectors
          [,1]       [,2]
[1,] 0.8889746 -0.7071068
[2,] 0.4579566  0.7071068

> y = x[,1]

这将产生以下错误消息：

Error in x[, 1] : incorrect number of dimensions

eigen 返回一个命名列表，其中特征值命名为值，特征向量命名为向量。访问列表的这个组件。我们使用美元符号。因此，这就是提取矩阵的代码 x$vectors 的原因。

【讨论】：