I am trying to use the functionality of Expectation
and Probability
for random graphs, in particular for percolation models.
For example, I would like to be able to compute the expected graph distance between two fixed vertices in a complete graph with independent exponential edge weights. The following naive code does not work.
n = 5;
Expectation[
GraphDistance[
CompleteGraph[n,EdgeWeight -> Table[ew[m], {m,1,(n (n-1))/2}]], 1, 2],
Table[ew[m] \[Distributed] ExponentialDistribution[1], {m,1,(n(n-1))/2}]
]//FullForm
yields
GraphDistance[Graph[List[1,2,3,4,5],List[Null,SparseArray[Automatic,List[5,5],0,List[1,List[List[0,4,8,12,16,20],List[List[2],List[3],List[4],List[5],List[1],List[3],List[4],List[5],List[1],List[2],List[4],List[5],List[1],List[2],List[3],List[5],List[1],List[2],List[3],List[4]]],Pattern]]],List[Rule[EdgeWeight,List[ew[1],ew[2],ew[3],ew[4],ew[5],ew[6],ew[7],ew[8],ew[9],ew[10]]],Rule[GraphLayout,"CircularEmbedding"]]],1,2]
indicating that the random weights are not actually assigned to the edges.
A similar example would be the computation of the probability that two fixed vertices in a complete graph are in the same connected component after edges are deleted with probability $1-p$ and retained with probability $p$. Again the following naive code does not work.
n = 5; p = 1/3;
Probability[
GraphDistance[
CompleteGraph[n, EdgeWeight -> Table[ew[m], {m, 1, (n(n-1))/2}]], 1, 2] == 0,
Table[ew[m] \[Distributed] EmpiricalDistribution[{1-p, p} -> {1, 0}], {m,1,(n(n-1))/2}]
]
Ideally, I would like to be able to compute this probability (a polynomial in $p$) for symbolic $p$. Of course I would also be grateful for other approaches to the problem, but using Expectation
and Probability
would feel very intuitive to me.
Update:
The problem can be circumvented by computing the probability of a generic event which is only specified after Probability
has been applied:
Clear[f, g]
n = 5;
f[ewarray_] :=
GraphDistance[CompleteGraph[n, EdgeWeight -> ewarray], 1, 2]
Probability[
g[Table[ew[m], {m, 1, (n(n-1))/2}]] == 0, Table[ew[m] \[Distributed] BernoulliDistribution[1-p], {m, 1, (n(n-1))/2}]
] /. g -> f // Simplify
The results are the same as @ubpdqn's MC calculations. I have experimented with Hold
and Release
to avoid the use of the extra function f
but do not know enough those functions to make it work.