Source code for calculation of stationary distribution in R

Question

Take a look at this link.

I am trying to understand the following source code meant for finding stationary distribution of a matrix:

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

I tested the following source code myself:

> 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
>

looks like the command

y <- x$vectors[,1]

is selecting the 1st column of the matrix.

Why wasn't that simply written like the following?

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

What was the reason for introduction of a dollar sign and vector keyword?

The dollar operator was inherited from S/Splus. Better avoid it. (your observation is correct, IMO) — wildplasser, May 10 '19 at 00:21
*"What was the reason for introduction of a dollar sign and vector keyword?"* I'm confused what you're asking. `eigen` returns a `list` with element `vectors` being the eigenvectors stored in a matrix. So you're simply accessing the `vectors` element from the `list` returned by `eigen`. There's nothing wrong with accessing `list` elements with `$` indexing, e.g. `lst <- list(a=1, b=2); lst$a`. — Maurits Evers, May 10 '19 at 00:21
@MauritsEvers, what do u think about **wildplasser**'s comment? — user366312, May 10 '19 at 00:24
@user366312 I think wildplasser may have misunderstood your question. As I said there's *absolutely* nothing wrong with using `$` to access elements from a named `list`. You'd use `$`-indexing to access elements (columns) from a `data.frame` as well, wouldn't you? — Maurits Evers, May 10 '19 at 00:26
@MauritsEvers, *I think wildplasser may have misunderstood your question.* -- in what way? I don't think so. — user366312, May 10 '19 at 00:41
Hmm. Perhaps I am misunderstanding the issue. Accessing named `list` elements using `$` notation is **absolutely standard practice in R**. For a specific example involving `eigen`, see e.g. [Eigenvalues and eigenvectors from the official "An Introduction to R" on CRAN](https://cran.r-project.org/doc/manuals/r-release/R-intro.html#Eigenvalues-and-eigenvectors); `$` indexing is certainly *not* something that should be avoided in R. — Maurits Evers, May 10 '19 at 01:34
I don't think your proposed simpler code produces the correct output. As mentioned, `eigen` returns a list with two objects. You have to tell R what item in the list you want (here `vectors`) THEN tell R what part of the matrix you want. `x[[2]][,1]` would be equivalent to `x$vectors[,1] `. Both of these are common. — Peter_Evan, May 10 '19 at 01:39
Seems very likely that @wildplasser is not a regular user of R. While it might be more safe to avoid `$colname` in functions, were one to do so, one would need to know that it is a shorthand for `[['colname']]` — IRTFM, May 10 '19 at 03:44

score 1 · Accepted Answer · answered May 11 '19 at 17:19

Let's test out your proposal:

> 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]

This would produce the following error message:

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

eigen returns a named list, with eigenvalues named values and eigenvectors named vectors. To access this component of the list. we use the dollar sign. Hence, that is why the code x$vectors which extract the matrix.

Source code for calculation of stationary distribution in R

1 Answers1