2

I am trying to generate a version of Watts' model of cascading shocks to a network following the discussion on p. 20-25 of this article.

The basic premise is to use the RandomGraph function to generate a graph, then create a table of characteristics on from VertexList. The values in the table (x) are the vertex number, the 'persuadability' of the vertex, and the current state (0,1). Everyone is at 0 to start, and we put in an initial shock of 5 vertices going from 0 to 1. Then if the ratio of neighbours who are 1s is above the 'persuadability' number, that neighbour changes its state from 0 to 1.

UPDATE 2

The code below does what I want, in that it iterates the cascade and records each step, and also does it for different initial system shocks.

The code does most of what I originally asked about -

Ultimately I would like to be able to do the following:

  1. Use the code to store a table of x tables at each iteration. This could be used to animate the graph later (using HighlightGraph[SubGraph[...]]).

  2. Use the code to simulate various random shocks on the same graph by iterating until the table x does not change.

  3. Maybe a nice Manipulate interface at some point.

However, the evaluation takes about 20 minutes, and it also iterates a fixed number of times iter instead of sufficient for convergence. So I have some specific questions.

  1. Why is the assignment record[[l]]=x necessary in the Table loop?
  2. How can I get the Table to iterate till record[[l]]==record[[l-1]]? (Mike suggests FixedPointList, but I can't seem to get the syntax to work)
  3. How can I write the Do loop with the If statement so I can use ParallelDo or similar? (much like the answer here, but with conditions from two lists)
  4. What are some other ways to speed this up generally?

UPDATE 3

This part of the code takes about half a second (which seems slow, since it is called every iteration, or 750 times)

nextn = Table[
          Drop[VertexList[NeighborhoodGraph[y, neighbourVL[[i]]]], 1], {i,
            1, Length[neighbourVL]}];

This generates a table with each row a list of neighbours of the vertexes in the list neighbourVL. Is there a way to get this list with a faster method?

Thanks everyone.

  (* Note 20min evaluation time *)

y = RandomGraph[{100, 150}];
edges = EdgeList[y];
nodes = Partition[
   Flatten[Riffle[VertexList[y], Table[{RandomReal[0.7], 0}, {100}]]],
    3];
iter = 15;
record = Table[0, {iter}];

output = Table[

   shock = Transpose[{RandomSample[Range[100], 5], Table[3, {5}]}];
   x = ReplacePart[nodes, shock -> 1];

   Table[
    converted = Cases[x, {_, _, 1}];

    neighbourVL = 
     Flatten[Table[
       Drop[VertexList[NeighborhoodGraph[y, converted[[i, 1]], 1]], 
        1], {i, 1, Length[converted]}]];

    (* Vertex list of neighbours' first neighbours to be evaluated
with If function *)

    nextn = Table[
      Drop[VertexList[NeighborhoodGraph[y, neighbourVL[[i]]]], 1], {i,
        1, Length[neighbourVL]}];

    (* Number of neighbours of vertex 1's first neighbour who are 1s *)

        neighbourcount = 
     Table[Table[
       Count[x, {Take[nextn[[i, j]]], _, 1}], {j, 1, 
        Length[nextn[[i, All]]]}], {i, 1, Length[nextn]}];

    threshold = 
     N[Total[neighbourcount, {2}]/
       Table[Length[nextn[[i, All]]], {i, 1, Length[nextn]}]];

    Do[

     If[threshold[[i]] > x[[neighbourVL[[i]], 2]], 
       x = ReplacePart[x, {neighbourVL[[i]], 3} -> 1]];, {i, 
      Length[nextn]}];


    record[[l]] = x

    , {l, 1, iter}]

   , {50}];


Histogram[Total[output[[All, -1, All, -1]], {2}]]

Manipulate[
 Table[HighlightGraph[y, 
   Subgraph[y, Cases[output[[j, i]], {_, _, 1}][[All, 1]]]]
  , {j, 1, 6}], {i, 1, iter, 1}]
Cameron Murray
  • 921
  • 7
  • 16
  • kguler, for me the {Manipulate} seems to work with or without the {With}. Excuse the beginner question, but what is the purpose of defining l as itself? – Cameron Murray Oct 23 '12 at 08:53
  • you are right, it does work with or without With. Sorry for the confusion. – kglr Oct 23 '12 at 10:01
  • The link doesn't go to an article by Watt. I found this link (http://oz.stern.nyu.edu/phd04/watts.pdf) to a PNAS article, but I'd like to see yours. Important topic. – George Wolfe Oct 23 '12 at 12:17
  • Thanks George. Your link is the original article. The link I have goes to a paper by Ormerod which simply had a very basic example of Watt's model. – Cameron Murray Oct 24 '12 at 00:11

1 Answers1

2

I have code that works about 50x faster now. The trick was to use the GraphUtilities package function NeighborhoodVertices[] rather than VertexList[NeighborhoodGraph[]]. It sped up these two lines from about 0.5sec each to 0.007sec. Since these are called 1000 times or more, that's a big help.

I am still interested in a more functional way of doing this. I have two Tables nested with each other, with the inner Table running a Do loop each iteration. Any pointers to previous questions about converted code to functional and faster code is appreciated (I think I have read more of them). It is probably the inner Table function that could be a separate function to map over a table of shocked node lists (x).

It takes about 30sec now (40 shocks with 20 iterations each). Would be great to improve that.

Here is the code

shockNum = 40;
iter = 20;
y = RandomGraph[{100, 130}];
nodes = Partition[
   Flatten[Riffle[VertexList[y], Table[{RandomReal[0.6], 0}, {100}]]],
    3];

output = Table[
   shock = Transpose[{RandomSample[Range[100], 5], Table[3, {5}]}];
   x = ReplacePart[nodes, shock -> 1];

   Table[
    converted = Cases[x, {_, _, 1}];
    neighbourVL = 
     Flatten[Table[
       Drop[NeighborhoodVertices[y, converted[[i, 1]], 1], 1], {i, 1, 
        Length[converted]}]];
    nextn = 
     Table[Drop[NeighborhoodVertices[y, neighbourVL[[i]], 1]], {i, 1, 
       Length[neighbourVL]}];
    neighbourcount = 
     Table[Table[
       Count[x, {Take[nextn[[i, j]]], _, 1}], {j, 1, 
        Length[nextn[[i, All]]]}], {i, 1, Length[nextn]}];
    threshold = 
     N[Total[neighbourcount, {2}]/
       Table[Length[nextn[[i, All]]], {i, 1, Length[nextn]}]];
    Do[
     If[threshold[[i]] > x[[neighbourVL[[i]], 2]], 
       x = ReplacePart[x, {neighbourVL[[i]], 3} -> 1]];, {i, 
      Length[nextn]}];
    x
    , {iter}]
   , {shockNum}];

Manipulate[
 Column[{
   Histogram[Total[output[[All, i, All, -1]], {2}], 10, 
    PlotRange -> {{0, 100}, {0, shockNum/2}}, ImageSize -> 550, 
    AspectRatio -> 0.3, AxesOrigin -> {0, 0}],
   Partition[
     Table[HighlightGraph[y, 
       Subgraph[y, Cases[output[[j, i]], {_, _, 1}][[All, 1]]]]
      , {j, 1, 6}], 3] // Grid}]
 , {{i, 1, "Iteration"}, 1, iter, 1, Appearance -> "Labeled", 
  ContinuousAction -> False}, 
 FrameLabel -> {{None, None}, {None, 
    "Watt's model of cascading network failure"}}, 
 AppearanceElements -> All, LabelStyle -> Directive[Black, Medium]]

And here is the output.

enter image description here

Cameron Murray
  • 921
  • 7
  • 16