Random Boolean Networks

Kai Willadsen and Ben Skellett

The aim of this tutorial is to introduce the structure and dynamics of Random Boolean Networks (RBNs). It will also provide an introduction to Matlab and the RBN toolbox for Matlab.

This document describes the mechanics and dynamics of Random Boolean Networks. Each section includes a series of questions. At the end of the document there is a series of tasks requiring the use of Matlab to investigate the structure and dynamics of RBNs.

For more information on RBNs see:

Motivation

The Random Boolean Network model was introduced by Kauffman in 1969 as an abstract model of the dynamics of gene regulation. In this model, nodes in the network correspond to genes, and edges between nodes indicate regulatory relationships between genes. That is, if there is an edge between node A and node B, then gene A regulates gene B. Attractors in the network (described below) have been argued to correspond to cell types, and the transient dynamics of the network have been likened to a developmental process. The aim of using Random Boolean Networks as a model of gene regulation is to ask questions about the nature of self-organisation in biological systems, and to explore the ways in which simple changes to a system's rules and biases can affect its global behavioural properties.

In addition to modelling the general properties of the dynamics of gene regulation, Random Boolean Networks are also prototypical discrete dynamical systems and are illustrative of combinatorial structure.

Mechanics

A Random Boolean Network consists of N randomly connected nodes, each of which has a binary state: on or off (1 or 0). In NK networks (a.k.a. Kauffman networks) every node receives exactly K inputs chosen randomly from other nodes in the network and in most cases self-connection is allowed (Figure 1 shows an example of the topology of an RBN).


Figure 1: Example of the topology of an RBN.


The state of each node in the network at time t+1 is determined by the states of its inputs at time t through a randomly generated Boolean function which can be represented with a look-up table for each node (see Table 1). Typically the Boolean functions do not change throughout the lifetime of the network. In an RBN the Boolean function for each node maps each of the 2K possible input combinations to an output state of 0 or 1.

Table 1: Example of a Boolean function with K inputs represented by a look-up table. Each node in an NK network has its own Boolean function with a randomly chosen output state for each of the 2K possible input combinations.


The network is given a random initial state by assigning each node a value of 0 or 1. The value of each node at the next time step is determined using the state of its inputs and each node's Boolean function. All nodes are updated at the same time (synchronously). The new state of the network is generated and then used to determine the next state and so on. An example of the network activity is displayed in Figure 2.

Figure 2: Activity pattern for an RBN with 16 nodes for 50 time steps. The initial state is the column furthest to the left with nodes represented vertically and time moving to the right.


Questions

  1. What are the values of N and K in Figure 1?
  2. Construct a random look-up table for a node with 3 inputs by arbitrarily choosing the outputs.
  3. What is the period of the expression pattern in Figure 2?

State Space

Since each node can be in one of two states, the total number of possible states for the entire network is 2N. Thus the system has a finite state space where each network state can be represented as a bit-string of length N. The use of deterministic Boolean functions means that the transitions between the network states are also deterministic and with a finite state space, the network will always fall into periodic behaviour.In other words, if the network is in state A and proceeds into state B at the next time step, then this transition will occur every time the system is in state A. As soon as one state is repeated, the network has entered a periodic attractor with a period length equal to the number of iterations since the network was last in same state. There may be a number of periodic attractors in the state space with period lengths ranging from 1 to a possible (though extremely unlikely) 2N.

Figure 3: A periodic attractor (yellow) and its basin of attraction (cyan). Each point in the state space represents a network state.


Not every state of the network is part of a periodic attractor. The other states lie in the basins of attraction for each attractor which consist of the transient states through which the network passes before it enters an attractor. Due to the deterministic dynamics none of the basins can overlap so the state space is divided into a number of periodic attractors with their associated basins of attraction.

Figure 4: The entire state space of an RBN with 10 nodes. Note: Self connections do not appear so a period-1 attractor appears to have no outputs although each network state must have exactly one output.


The number of attractors and the size of their periods are important quantities for RBNs. If a network is in a period 1 attractor its behaviour will be fixed. Behaviour with a low period is also easily recognisable as periodic (see Figure 2), however, if the behaviour has a very high period, the activity pattern will appear random as it may take many thousands of time steps before the pattern repeats, even for modest sized networks. Long periods typically occur in networks that can be classified as chaotic (see below) and vice versa.

Questions

  1. If an RBN consists of 30 nodes, how many distinct network activity patterns are there (in decimal)? What about an RBN with 16 nodes or 10 nodes?
  2. What is the period of the attractor in Figure 3?
  3. How many attractors are there in Figure 4?
  4. What is the average period length of the attractors in Figure 4?
    Hint: Follow the arrows. Self connections do not appear so a period-1 attractor appears to have no outputs. Period-2 attractors appear to have bi-directional arrows: ↔

Network Classification

Network behaviour can be classified into 3 regimes: frozen, chaotic and critical. The existence of chaotic behaviour is identified by the way in which perturbations (e.g. changing the state of a node) spread through the network. To measure such damage spreading the activity pattern of an RBN is recorded, this is the original activity pattern. Then the initial network state is perturbed by changing the state of one node from a 0 to a 1 or vice versa. The perturbed activity pattern is monitored over time. After a number of iterations, if the damage has spread so that the network state now differs from the original activity pattern in an increasing number of places, the network is called chaotic — characterized by sensitive dependence on initial conditions. If the perturbation has disappeared so that both networks display the same pattern, the network is called frozen (or ordered). Alternatively, the states of the two networks still differ in a limited number of places so that the damage may not have spread or disappeared, in this situation the network is said to display critical behaviour.

The key factor in determining network behaviour is the connectivity. Averaged over many networks, critical behaviour occurs in NK networks for K = 2, frozen behaviour occurs for K < 2 and chaotic behaviour occurs for K > 2. It is important to understand that this classification applies when averaged over many networks and that an individual network with K > 2 can display ordered dynamics with low period attractors.

An associated property in the classification of networks is the way in which period length grows with increasing the number of nodes (see Table 2).

Regime K Period Length Damage Spreading Matter Analogy
Frozen

Less than 2

Independent of N (typically period 1)

Perturbations die out; stable dynamics.

Solid
Critical

2

Grows algebraic with N (believed to be proportional to N1/2).

Perturbations continue to affect a limited number of nodes; marginally stable dynamics.

Transitional (Liquid?)
Chaotic

Greater than 2

 Grows exponentially with N. Perturbations cause extensive changes to the network state; unstable dynamics. Gas

Table 2: RBN Classification


A list of the properties of Random Boolean Networks outlined above is available in the Summary of Random Boolean Network Properties.

Tasks

Setup

The following tasks are designed to be done using the Matlab mathematical programming environment. If you are not familiar with Matlab, you should first look at the One Page Introduction to MATLAB. You will also be using a set of tools called the Matlab Random Boolean Network Toolbox.

  1. Download the RBN Toolbox from http://fias.uni-frankfurt.de/~willadsen/RBN/RBNToolbox.zip
  2. Unzip the file
    Just click the 'Next' button until the wizard finishes. Extract all of the files to your "My Documents" directory
  3. Start Matlab
  4. Change to the directory in which you unzipped the Toolbox
    Click the '...' button next to the 'Current directory' entry in the toolbar. Go to "C:\Documents and Settings\accs<number>\My Documents"

Documentation for Matlab (and RBN Toolbox) commands can be obtained by typing help <command> at the Matlab prompt, where <command> is the function that you wish to access help for. Documentation is also available at the toolbox website.

Create and investigate a Random Booleam Network

  1. Create a Random Boolean Network and view the network's structure
    1. Specify N (number of nodes in the network) and K (number of inputs per node)
      N = 10; K = 2;
    2. Set the seed of the network
      rand ('state', 12345);
    3. Create a Random Boolean Network
      [node, conn, rules] = bsn (N, K, 'arrow');
      The node variable contains the set of nodes in the network, the rules variable contains the boolean updating functions for the network. The conn variable contains the connectivity or adjacency matrix of the network. A network with N nodes can be represented by an N by N adjacency matrix. The matrix, conn, is a matrix of zeros and ones in which, if there is an edge from node i to node j, conn(i,j) is 1, and 0 otherwise.
  2. Investigate the in-degree and out-degree of the network
    1. Create a histogram of the in-degree distribution of the network
      hist (sum (conn), 1:max(sum(conn)));
    2. Create a histogram of the out-degree distribution of the network
      hist (sum (conn'), 1:max(sum(conn')));
  3. Run the network to obtain an expression pattern and find the period of an attractor
    1. Run the network and display the expression pattern
      [node, tsm] = displayEvolution (node, 50, 'CRBN');
    2. Count the number of states in the attractor that the network reaches. To check your answer run,
      [attrLength, attrStates] = findAttractor (tsm);
      attrLength
  4. Generate the state space for the network and export it to Pajek for visualisation
    1. To generate the state space for your RBN run,
      stateSpace = fillStateSpace (node);
      Warning: This enumerates every possible state in the network, taking 2N space and time. If you increase N beyond 10, expect to wait a long time.
    2. Export the state space adjacency matrix to Pajek using,
      writeStateSpaceToPajek (stateSpace, 'myRBN');
    3. Start Pajek and open 'myRBN.net' from the directory in which Matlab was working.
    4. Use draw-partition to visualise the network and change the layout using Layout → Energy → Kamada-Kawai → Free.
    5. In Pajek, compute the number of weakly-connected components. Use Draw → partition to visualise the network. Compute the Partition → Core → Input to highlight the attractors.

Explore the different classes of network behaviour

  1. Find a Random Boolean Network that demonstrates ordered behaviour.
    1. Altering the network connectivity K and the random number seed, find a network that shows a point-attractor.
      Note that when you change the random number seed, you need to recreate the network using the bsn command as above.
    2. Check that the network displays a point attractor for several different starting states. You can set the starting state of the network using,
      node = setStateVector(node, stateVector)
      where stateVector is a vector of Booleans representing whether each node in the network is on or off. If you want to find out what the current state of the network is, you can use,
      state = getStateVector(node)
      You can then modify this state by setting individual elements,
      state(4) = 1
      and then set the network activation to the new state, node = setStateVector(node,state)
    3. Export the network's state space to Pajek. Make sure that you give the Pajek file a meaningful name, because you'll be comparing different state spaces later.
  2. Find a Random Boolean Network that demonstrates cyclic behaviour.
    1. Using the same procedure as above, find a network that consistently demonstrates cyclic behaviour and export its state space to Pajek
  3. Find a Random Boolean Network that demonstrates chaotic behaviour.
    1. Using the same procedure as above, find a network that consistently demonstrates chaotic behaviour and export its state space to Pajek
  4. Compare ordered, cyclic and chaotic state spaces
    1. One at a time, load the saved state spaces into Pajek.
    2. Visualise the state spaces as above.
    3. Identify visual differences between them (other than the length of the attractors). Can you infer anything about network behaviour from these differences?