spa.py

""" Generate graphs from the spatial preferential attachment model,
and animate them"""

from pymc import runiform, rnormal
from networkx import DiGraph, draw
from pylab import clf, axis, savefig
from numpy import maximum, minimum


def perturb(x, e):
    return maximum(0., minimum(rnormal(x, e**-2), 1.))


def distance(x, y):
    return max(abs(x-y))


def make_graph(n, d, p, A1, A2, A3, position=None, affinity=None):
    """ Create a directed graph according to the SPA model, using the
    randomness provided in position and affinity, or generating it if
    necessary.
    """
    if position == None:
        position = runiform(0, 1, size=(n, d))
    if affinity == None:
        affinity = runiform(0, 1, size=(n,n))
        
    G = DiGraph()
    G.add_nodes_from(range(n))

    for t in range(1, n):
        for s in range(t):
            if distance(position[s], position[t]) \
                    <= (A1 * G.in_degree(s) + A2) / (t + A3):
                if affinity[t, s] <= p:
                    G.add_edge(t, s)

    return G



class SPA:
    def __init__(self, n=100, d=2, p=.3, A=[.1,5.,.1], delta=.01):
        # set dimensions
        self.n = n
        self.d = d

        # set model parameters
        self.p = p
        self.A1, self.A2, self.A3 = A

        # set time stuff
        self.t = 0
        self.delta = delta
        
        # initialize parameters
        self.position = runiform(0, 1, size=(self.n, self.d))
        self.affinity = runiform(0, 1, size=(self.n, self.n))

    def update(self):
        self.position = perturb(self.position, self.delta)
        self.affinity = perturb(self.affinity, self.delta)
        
        self.G = make_graph(self.n, self.d, self.p,
                            self.A1, self.A2, self.A3,
                            position=self.position,
                            affinity=self.affinity)
        clf()
        draw(self.G, pos=self.position, alpha=.8)
        axis([0,1,0,1])

        savefig('spa_%03d.png' % self.t)

        self.t += 1

if __name__ == '__main__':
    s = SPA()
    for i in range(900):
        s.update()
    import subprocess
    subprocess.call('mencoder mf://*.png -mf w=800:h=600 -ovc x264 -of avi -o spa.avi',
                    shell=True)