The `rviewgraph`

package provides functions for creating,
changing, and viewing a graph. It originally contained only one function
called `rViewGraph()`

that allowed viewing a fixed graph
using an animated graphical user interface (GUI). A newer function
called ’vg()` is now available and allows for viewing a graph the
structure of which can change dynamically. To illustrate this we will
implement a simple Markov chain Monte Carlo (MCMC) sampler with a given
target distribution.

First ensure that `rviewgraph`

is available, and attempt
to start both an old version of the GUI using `rViewGraph()`

,
and new version of the GUI using `vg()`

. For the old version
we need to specify a graph in advance so we generate some random edge
end points.

```
library(rviewgraph)
#> Loading required package: rJava
= rViewGraph(sample(1:10,10,TRUE),sample(1:10,10,TRUE))
oldgui #> [1] "rViewGraph only runs in interactive mode."
is.null(oldgui)
#> [1] TRUE
if (interactive() && !is.null(oldgui))
{oldgui.stop()
oldgui.hide()
}
= vg()
v #> [1] "Running without a visible GUI."
is.null(v)
#> [1] FALSE
```

We see that `rViewGraph()`

returns `NULL`

and
prints an error message. This is because the old programs can only be
run in interactive mode. The new version, `vg()`

, however,
will return a list of functions that can be run non-interactively, but,
of course, without a visible GUI. As this vignette was produced
non-interactively, there is no GUI, but if you are replicating it from a
live `R`

session you should see a new `Java`

window that has appeared with a blank canvas on which the graph will be
viewed. There will also be a GUI showing the random graph generated in
the call to `rViewGraph()`

, so we remove that.

The returned value `v`

is a list of functions that we will
use to display the graph as we run the sampler.

First we define a probability mass function on graphs, in fact we need only specify this up to a constant of proportionality. This simple example depends only on the number of edges in the graph, and as such it would be more efficient to track this count as we iterate through the MCMC scheme, however, this will illustrate how the structure of the graph can be queried, and also makes the MCMC scheme more easily adapted to different probability distributions. For numerical stability, we work with logs of probabilities.

This probability distribution favours more sparsely edged graphs if
we set `a > 0`

, and more densely edged ones if
`a < 0`

. We can query the graph for a two dimensional
array that specifies the edges. Because this graph is undirected, each
edge is returned twice, once in each direction, so we divide to get the
number of edges.

```
= function(g)
edges
{dim(g$getEdges())[1] / 2
}
= function(g,a)
logProb
{-edges(g)*a
}
```

Here is a function that flips an edge in a graph. That is, if a pair
of vertices are connected it disconnects them, if disconnected, it
connects them. This will be our Metropolis proposal scheme, and it is
reversible. The argument `g`

is assumed to be a list that was
returned by calling `vg()`

.

```
= function(g,x,y)
flip
{if (g$connects(x,y))
$disconnect(x,y)
gelse
$connect(x,y)
g }
```

Now we define a Metropolis MCMC proposal and acceptance step. The
calls to `Sys.sleep()`

just slow the process down to make the
process easier to follow in the viewer when used interactively. The
colour changes will show the proposals being made.

```
= function(g, a, t = 1)
step
{# Find the probability of the current graph.
= logProb(g,a)
delta
# Choose a pair of random vertices to flip, and display them.
= sample(g$getVertices(),2,replace=FALSE)
x $colour(x,"red")
gSys.sleep(t/3)
# Then make the flip.
flip(g,x[1],x[2])
Sys.sleep(t/3)
# Find the ratio of the probabilities of the proposed new graph
# and the previous incumbent.
= logProb(g,a) - delta
delta
# Then with probability equal to the minimum of 1 and this ratio
# accept the new graph by doing nothing, otherwise reject it by
# undoing the flip.
if (log(runif(1)) > delta)
flip(g,x[1],x[2])
# Reset the colours, and pause before the next step.
$colour(x)
gSys.sleep(t/3)
}
```

We can now combine these functions to make a complete MCMC scheme
using here a graph of `n = 20`

vertices. As an example of an
output, we can calculate the mean number of edges over the simulation.
If we set the penalty to 0, we get a random walk over graphs of size
`n`

, so the simulated mean number of edges should be about
`n(n-1)/4`

.

If this is run in interactive mode, this is primarily to illustrate how the GUI shows the progress of the sampler, so we set the sleep parameter to give one step per second and make 50 updates. If run non-interactively, we don’t need to slow down the sampler to see the progress, so we make more iterations and include a burn in period to allow the sampler to move away from the initial state and obtain better estimates of the mean number of edges.

The scheme is run using the viewer that was created at the beginning of this vignette.

```
# To make the simulation reproducible.
set.seed(1)
= 20
n $add(1:n)
v
= 0
pen
if (interactive())
{= 1
wait = 50
s else
}
{= 0
wait = 1000
s for (i in 1:s)
step(v,pen,wait)
}
= 0
medge for (i in 1:s)
{step(v,pen,wait)
= medge + edges(v)
medge
}c(medge/s, n*(n-1)/4)
#> [1] 98.164 95.000
```

A positive penalty should decrease the mean number of edges.

```
= 1
pen = 0
medge for (i in 1:s)
{step(v,pen,wait)
= medge + edges(v)
medge
}/s
medge#> [1] 58.688
```

And a negative one should increase the mean number of edges.

```
= -1
pen = 0
medge for (i in 1:s)
{step(v,pen,wait)
= medge + edges(v)
medge
}/s
medge#> [1] 127.254
```