Asymptotic Mutual Information 
between Boolean Net Regions

Tommaso Bolognesi
CNR/ISTI, Pisa

ABSTRACT

Boolean networks (‘bool nets’) are sequential, dynamical systems introduced by S. Kauffman for modelling biological, genetic systems.  They are a generalisation of cellular automata.  
In this talk we conceive bool nets as composed of two sub-nets P and Q that interact by sharing  state variables, with the degree of coupling controlled by the user.
Some interactive Mathematica demonstrations are introduced.  Their purposes are: 
(i) to specify and simulate synchronous, simple or composite bool nets; 
(ii) to collect statistical data on the mutual information M(P,Q) between the two bool net components P and Q;
(iii) to relate M(P,Q) to the degree of coupling between P and Q, and to the computation length.

Context and motivations

Seminal work by Tononi, Balduzzi, Oizumi, Albantakis (since 2004) on Integrated Information (Φ) Theory  ⟺ consciousness

Φ applies to discrete dynamical systems -- typically boolean nets -- and quantifies  ‘the whole (when parts interact) is more than the sum of the parts (no interaction)’

How does M(P, Q) -- Mutual Information between system parts  -- fit into the Φ picture?

Code

randJoint2d[n_,m_]:=With[{nums=RandomInteger[{1,100},n*m]},Partition[nums/Total[nums],m]]

randomBoolNet[n_,k_]:=
spec[
(*boolFunCodeList*)
RandomInteger[{0,2^(2^k)-1},{n}],
(*argIndexTupleList*)
Table[RandomSample[Range[n],k],{n}]]

coupInducer[nn_,k_,argIndexTupleList_,coup_]:=
With[{positions=RandomSample[Range[1,nn*k],coup]},
Partition[
ReplacePart[
Flatten[argIndexTupleList],
Table[p->If[p>nn*k/2,RandomInteger[{1,nn/2}],RandomInteger[{nn/2+1,nn}]],{p,positions}]],
k]]

randomBoolNetOfCoupMULTIEDGE[nn_,k_,coup_]:=
With[{bipartRbn=spec[
(*boolFunCodeList*)
RandomInteger[{0,2^(2^k)-1},{nn}],
(*argIndexTupleList*)
Join[
Table[RandomChoice[Range[1,nn/2],{k}],{nn/2}],
Table[RandomChoice[Range[nn/2+1,nn],{k}],{nn/2}]]]},
ReplacePart[bipartRbn,2->coupInducer[nn,k,bipartRbn[[2]],coup]]]

boolNetGraph[bnSpec_]:=With[{fromNodesXnodePairs=Thread[List[bnSpec[[2]],Range[1,Length[bnSpec[[2]]]]]]},
Graph[Range[1,Length[bnSpec[[1]]]],Flatten[Thread[DirectedEdge[#[[1]],#[[2]]]]&/@fromNodesXnodePairs](*,VertexLabels->"Name"*)]]

nodeCoordInterpol[{coordInit_,coordFin_},i_]:=(* i = max 40th step *)
{coordInit[[1]]+i*(coordFin[[1]]-coordInit[[1]])/50,
coordInit[[2]]+i*(coordFin[[2]]-coordInit[[2]])/50}

nodesCoordInterpol[coordsInit_,coordsFin_,i_]:=
nodeCoordInterpol[#,i]&/@Transpose[{coordsInit,coordsFin}]

(* code in [0..2^(2^nArgs)-1] (function code) *)
boolFunNargs[code_,bitNtuple_]:=IntegerDigits[code,2,2^Length[bitNtuple]][[FromDigits[bitNtuple,2]+1]]

boolNetNargsStep[bitTupleState_,boolFunCodeList_,argIndexTupleList_]:=
Table[boolFunNargs[boolFunCodeList[[i]],bitTupleState[[argIndexTupleList[[i]]]]],{i,1,Length[bitTupleState]}]

boolNetNargsStateGraph[boolFunCodeList_,argIndexTupleList_]:=
Graph[
Prepend[Table[i,{i,1,2^Length[boolFunCodeList]-1}],
Labeled[0,Row[{
(*000*)Style[IntegerString[0,2,Length[boolFunCodeList]/2],Blue],
(*000*)Style[IntegerString[0,2,Length[boolFunCodeList]/2],Green]}]]],
Table[DirectedEdge[i,FromDigits[boolNetNargsStep[IntegerDigits[i,2,Length[boolFunCodeList]],boolFunCodeList,argIndexTupleList],2]],{i,0,2^Length[boolFunCodeList]-1}],
ImageSize->220(*,VertexLabels->"Name"*)]

boolNetNargsStateGraph[boolNetSpec_]:=boolNetNargsStateGraph[boolNetSpec[[1]],boolNetSpec[[2]]]

(* state graph nodes: 0 to 2^n -1; TPM rows/cols: ordered from 1 to 2^n; 
NO ROWS OF ZERO'S s for flat Bool Nets) *)
fromStateGraphToTPM[gg_]:=
With[{n=Max[VertexList[gg]]+1},
If[Total[#]==0, (* # is a TPM row *)
#,
#/Total[#]]&/@
Fold[ReplacePart[#1,{#2[[1]]+1,#2[[2]]+1}->#1[[#2[[1]]+1,#2[[2]]+1]]+1]&,Table[0,{n},{n}],EdgeList[gg]]]

yPQstarApprox[nn_,triples_]:=
Fold[ReplacePart[#1,#2[[1]]->(#2[[3]]/2^nn)(1/#2[[2]])]&,Table[0,{2^nn}],triples]

triplesForCycle[cycle_,wccs_]:=
With[{cycleNodes=Union[Flatten[List@@@cycle]]},
Table[{n+1,Length[cycle],Length[Select[wccs,MemberQ[#,n]&][[1]]]},{n,cycleNodes}]]

triplesForSinks[sinks_,wccs_]:=
Table[{s+1,1,Length[Select[wccs,MemberQ[#,s]&][[1]]]},{s,sinks}]

triplesForCyclesAndSinks[cycles_,sinks_,wccs_]:=
Join[Flatten[Table[triplesForCycle[cycle,wccs],{cycle,cycles}],1],triplesForSinks[sinks,wccs]]

marginal[pxy_,xORy_]:=
Switch[xORy,
"x",Total/@pxy,
"y",Total/@Transpose[pxy]]

mutInf[pxy_]:=
With[{px=marginal[pxy,"x"],py=marginal[pxy,"y"]},
Total[
Flatten[
Table[If[pxy[[i,j]]==0,0,pxy[[i,j]]Log[2,pxy[[i,j]]/(px[[i]]*py[[j]])]],{i,1,Length[px]},{j,1,Length[py]}]]]]

syncMutInfoApprox[stateGG_]:=
With[{
WCCs=WeaklyConnectedComponents[stateGG],
CCs=ConnectedComponents[stateGG]},
With[{
nn=Log[2,VertexCount[stateGG]],
cycles=FindCycle[stateGG,Max[Length/@CCs],Length[WCCs]],
sinks=Union[Flatten[List@@@Cases[EdgeList[stateGG],x_\[DirectedEdge]x_]]]},
N[mutInf[Partition[yPQstarApprox[nn,triplesForCyclesAndSinks[cycles,sinks,WCCs]],2^(nn/2)]]]]]

1. Random Boolean Networks (RBN)

Sequential dynamical system - generalization of Cellular Automata

Introduced by S. Kauffman for modeling genetic regulatory networks [see Refs].

Manipulate[
SeedRandom[sr];
With[{
rbn=randomBoolNet[nn,k],
initBitTuple=RandomInteger[1,nn]},
With[{
rbnGraph=boolNetGraph[rbn],
bitTupleList=NestList[boolNetNargsStep[#,rbn[[1]],rbn[[2]]]&,initBitTuple,6]},
With[{
coordsInit=Table[PropertyValue[{rbnGraph,n},VertexCoordinates],{n,1,nn}],
coordsFin=Join[Table[{0,y},{y,0,1,2/(nn-1)}],Table[{1,y},{y,0,1,2/(nn-1)}]]},
Pane[
Graph[
VertexList[rbnGraph],
EdgeList[rbnGraph],
VertexLabels->
Table[i->Row[{i,", ", 
Tooltip[
Style[Row[{"fun ",rbn[[1,i]]}],Red],
Column[Thread[Rule[Table[IntegerDigits[j,2,k],{j,0,2^k-1}],IntegerDigits[rbn[[1,i]],2,2^k]]]]] }],{i,nn}],
VertexSize->.16,VertexStyle->Table[i->(bitTupleList[[step+1,i]]/.{0->White,1->Black}),{i,nn}],
ImageSize->430,
PlotLabel->Style[Row[{Superscript[y,step]," = ",bitTupleList[[step+1]]}],Large,Background->Lighter[Yellow]]],
ImageSize->{550,450},Scrollbars->True,Alignment->Center]]
]],
{{nn,4,"Num. of bool net nodes (nn)"},{3,4,5}},
{{k,2,"Arity of bool functions (k)"},{2,3},ControlType->SetterBar},
{{sr,10006,"random seed"},10000,20000,1},
{{step,0},0,6,1,ControlType->SetterBar}]

2. Bool Net bipartition

Split node set into subsets P and Q

P and Q interact by shared nodes

Coupling factor

(Subscript[y, P]^h) and (Subscript[y, Q]^h) as random variables

Manipulate[
SeedRandom[sr];
If[coup>nn*k,coup=nn*k];
With[{
rbn=randomBoolNetOfCoupMULTIEDGE[nn,k,coup],
initBitTuple=RandomInteger[1,nn]},
With[{
rbnGraph=boolNetGraph[rbn],
bitTupleList=NestList[boolNetNargsStep[#,rbn[[1]],rbn[[2]]]&,initBitTuple,6]},
With[{
coordsInit=Table[PropertyValue[{rbnGraph,n},VertexCoordinates],{n,1,nn}],
coordsFin=Join[Table[{0,y},{y,0,1,2/(nn-1)}],Table[{1,y},{y,0,1,2/(nn-1)}]]},
Pane[
Row[{If[unravel>42,Style["P",Bold,Blue,Large],Null],
Dynamic[
Graph[
VertexList[rbnGraph],
EdgeList[rbnGraph],
VertexCoordinates->nodesCoordInterpol[coordsInit,coordsFin,unravel],
VertexLabels->"Name",VertexSize->.16,VertexStyle->Table[i->(bitTupleList[[step+1,i]]/.{0->White,1->Black}),{i,1,nn}],
ImageSize->430,
PlotLabel->If[unravel>42,
Row[{
Style[Row[{Superscript[Subscript[y,P],step]," = ",Take[bitTupleList[[step+1]],nn/2]}],Large,Blue,Background->Lighter[Yellow]],
"       ",
Style[Row[{Superscript[Subscript[y,Q],step]," = ",Take[bitTupleList[[step+1]],-nn/2]}],Large,Green,Background->Lighter[Yellow]]}],
Null]]],
If[unravel>42,Style["Q",Bold,Green,Large],Null]},ImageMargins->0],
ImageSize->{550,600},Scrollbars->True,Alignment->Center]]
]],
{{nn,6,"Num. of bool net nodes (nn)"},{4,6,8,10}},
{{k,2,"Arity of bool functions (k)"},2,5,1,ControlType->SetterBar},
{{unravel,0,"P-Q unravel"},0,50,1},
{{coup,5,"P-Q Coupling factor (coup)"},0,nn*k,1},
{{sr,12345,"random seed"},10000,20000,1},
{{step,0},0,6,1,ControlType->SetterBar}]

3. Mutual information M(x,y) between generic random vars x and y

MutInfExpr=TraditionalForm[Sum[(Subscript[P, xy])[(Subscript[x, i]),(Subscript[y, j])](Subscript[Log, 2])[(Subscript[P, xy])[(Subscript[x, i]),(Subscript[y, j])]/((Subscript[P, x])[(Subscript[x, i])](Subscript[P, y])[(Subscript[y, j])])],{i,1,n},{j,1,m}]];
Manipulate[
SeedRandom[randSeed];
With[{joint=randJoint2d[n,m]},
Grid[
{{Grid[{{Style[Subscript["Joint P",xy],Large],Style[y,Large]},{Style[x,Large],MatrixForm[joint]}}],
Column[{
Style[Subscript["Marginal P",x],Large],
Column[Total/@joint,Spacings->1]},Spacings->2]},
{Column[{
Style[Subscript["Marginal P",y],Large],
Row[Total/@Transpose[joint],"   "]},Spacings->2],
Column[{
Style["Mutual Info M(x,y)=",Large],
Style[Row[{MutInfExpr, "="}],Blue,Bold],
Style[N@mutInf[joint],Blue,Large]
},Background->Lighter[Yellow], Spacings->2]}},
Frame->All,ItemSize->{{32,24},{30,20}},Alignment->{Center,Center}
]],
{{randSeed,1234,"random seed"},1000,2000,1},
{{n,8},2,8,1,ControlType->SetterBar},
{{m,8},2,8,1,ControlType->SetterBar}]

Venn diagram: (conditional) entropies H and mutual info. M

M(x,y) as M((Subscript[P, xy])) --- a function of a joint distribution

4. Sequence of state distributions

Transition Probability Matrix (TPM) for bool net PQ

SeedRandom[9094796]
nn=4;
rbn=randomBoolNet[nn,3];
net=boolNetGraph[rbn];
gg=boolNetNargsStateGraph[rbn];
tpm=fromStateGraphToTPM[gg];
Grid[{
{"Bool net PQ","Global state graph of PQ","TPM"},
{net,gg,ArrayPlot[tpm]}
},Frame->All]

(Subscript[Y, PQ]) = (Subscript[Y, PQ]^0)->  (Subscript[Y, PQ]^1)-> ...   (Subscript[Y, PQ]^h)->  ...

(Subscript[Y, PQ]^0) = {(1/16),(1/16),(1/16),(1/16),(1/16),(1/16),(1/16),(1/16),(1/16),(1/16),(1/16),(1/16),(1/16),(1/16),(1/16),(1/16)}
(Subscript[Y, PQ]^(h+1))= (Subscript[Y, PQ]^h).TPM

(Subscript[Y0, PQ])=Table[1/(2^nn),{2^nn}];
(Subscript[Y, PQ])=NestList[#.tpm&,(Subscript[Y0, PQ]),30];
ArrayPlot[(Subscript[Y, PQ])]

5. Sequence of joint distributions: (Subscript[Y, P])(Subscript[Y, Q]) = (Subscript[Y, P]^0)(Subscript[Y, Q]^0) -> (Subscript[Y, P]^1)(Subscript[Y, Q]^1) -> ... (Subscript[Y, P]^h)(Subscript[Y, Q]^h) -> ...
   and associated M values

from distrib. (Subscript[Y, PQ]^h) to joint distrib. (Subscript[Y, P]^h)(Subscript[Y, Q]^h)

YpYq=Partition[#,2^(nn/2)]&/@(Subscript[Y, PQ]);
Column[{
Row[Table[Column[Append[IntegerDigits[i,2,nn],"↓"]],{i,0,2^nn-1}]," "],
ArrayPlot[{(Subscript[Y, PQ])[[-1]]},ImageSize->450,Mesh->All]},
Center]

Grid[{
{"",Row[Table[Column[IntegerDigits[i,2,nn/2]],{i,0,2^(nn/2)-1}],"  "]},
{Column[Table[Row[IntegerDigits[i,2,nn/2]],{i,0,2^(nn/2)-1}],Spacings->1],ArrayPlot[YpYq[[-1]],ImageSize->100]}
}]

sequence of M((Subscript[Y, P]^h)(Subscript[Y, Q]^h)) values

ListPlot[mutInf/@YpYq,Joined->True,PlotRange->nn/2,AxesLabel->{Style["h",Large,Italic],Style["M",Large,Italic]}]

6. Approximating M((Subscript[Y, P]^h)(Subscript[Y, Q]^h)) using the attractors of the global graph

Manipulate[
SeedRandom[sr];
If[k>nn,k=nn];
With[{rbn=randomBoolNet[nn,k]},
With[{stateGG=boolNetNargsStateGraph[rbn]},
Grid[{{
Row[{"Global state graph and ",Text[Style["attractors",Red,FontFamily->"Courier"]]}],
Row[{"M values and ",Text[Style["approximation",Red,FontFamily->"Courier"]]}]},
{Pane[
HighlightGraph[
stateGG,
Join[
(*CYCLES*) Flatten[FindCycle[stateGG,Max[Length/@ConnectedComponents[stateGG]],Length[WeaklyConnectedComponents[stateGG]]],1],
(*SINKS*) Union[Flatten[List@@@Cases[EdgeList[stateGG],x_\[DirectedEdge]x_]]]]]
,ImageSize->{240,260},Scrollbars->True,Alignment->Center],
With[{tpm=fromStateGraphToTPM[stateGG]},
ListPlot[mutInf/@(Partition[#,2^(nn/2)]&/@NestList[#.tpm&,Table[1/(2^nn),{2^nn}],steps]),Joined->True,PlotRange->nn/2,AxesLabel->{Style["h",Large,Italic],Style["M",Large,Italic]},
GridLinesStyle->Directive[Red,Dashed,Thick],GridLines->{{},{syncMutInfoApprox[stateGG]}}]]}},
Frame->All,ItemSize->{{20,20},{4,29}},Alignment->{Center,Center},Spacings->{0,0}]]],
{{nn,4,"Num. of bool net nodes (nn)"},{4,6,8}},
{{k,3,"Arity of bool functions (k)"},2,8,1,ControlType->SetterBar},
{{steps,30},10,60,1},
{{sr,9094796,"random seed"},9094000,9095000,1},
SaveDefinitions->True]

The trick

as h ⟶ ∞,  distribution (Subscript[Y, PQ]^h)’selects’ the GG states in the attractors (Subscript[A, 1])... (Subscript[A, n])

One-to-one corresp. (Subscript[A, i]) ⟷ (Subscript[WCC, i]) (Weakly Conn. Components)

build (Subscript[Y, PQ]^( * )) such that:

Weight((Subscript[A, i])) ∼ Size((Subscript[WCC, i]))

all states in attractor (Subscript[A, i]) as equally likely (i = 1...n)

(Subscript[Y, PQ]^( * )) ⟶ (Subscript[Y, P]^( * ))(Subscript[Y, Q]^( * )) ⟶ M((Subscript[Y, P]^( * ))(Subscript[Y, Q]^( * )))

7. [Wolfram Demo] Asymptotic Mutual Information (Subscript[M, coup])((Subscript[y, P]^( * )), (Subscript[y, Q]^( * )))

as number of steps h ---> ∞

as a function of coupling factor ‘coup’

Manipulate[
SeedRandom[sr];
If[coup>nn*k,coup=nn*k];
With[{rbn=randomBoolNetOfCoupMULTIEDGE[nn,k,coup]},
With[{rbnGraph=boolNetGraph[rbn]},
With[{stateGG=boolNetNargsStateGraph[rbn]},
With[{
coordsInit=Table[PropertyValue[{rbnGraph,n},VertexCoordinates],{n,1,nn}],
coordsFin=Join[Table[{0,y},{y,0,1,2/(nn-1)}],Table[{1,y},{y,0,1,2/(nn-1)}]]},
Grid[{{
"Random bool net",
Row[{"Global state graph and ",Text[Style["attractors",Red,FontFamily->"Courier"]]}]},
{Row[{If[unravel>42,Style["P",Bold,Blue],Null],
Dynamic[
Graph[
ReplacePart[VertexList[rbnGraph],{1->Labeled[1,"1"],nn->Labeled[nn,nn]}],
EdgeList[rbnGraph],
VertexCoordinates->nodesCoordInterpol[coordsInit,coordsFin,unravel],(*VertexLabels->"Name",*)ImageSize->180]],
If[unravel>42,Style["Q",Bold,Green],Null]},ImageMargins->0],
Pane[
HighlightGraph[
stateGG,
Join[
(*CYCLES*) Flatten[FindCycle[stateGG,Max[Length/@ConnectedComponents[stateGG]],Length[WeaklyConnectedComponents[stateGG]]],1],
(*SINKS*) Union[Flatten[List@@@Cases[EdgeList[stateGG],x_\[DirectedEdge]x_]]]]]
,ImageSize->{240,260},Scrollbars->True,Alignment->Center]},
{Row[{Style["P",Bold,Blue],"-",Style["Q ",Bold,Green],"coupling factor 'coup' = ",coup}],Item[Row[{"MutualInfo(",Style["P,",Bold,Blue],Style["Q",Bold,Green],") = ",syncMutInfoApprox[stateGG]}],Background->Lighter[Yellow, 0.5]]}},
Frame->All,ItemSize->{{20,20},{4,29,4}},Alignment->{Center,Center},Spacings->{0,0}]]]
]],
Row[{
Control[{{nn,6,"Num. of bool net nodes (nn)"},{4,6,8,10}}],
Control[{{k,3,"Arity of bool functions (k)"},2,10,1,ControlType->SetterBar}],
Control[{{coup,5,"P-Q Coupling factor (coup)"},0,nn*k,1}],
Control[{{unravel,0,"P-Q unravel"},0,50,1}],
Control[{{sr,12345,"random seed"},10000,20000,1}]},"      "],
SaveDefinitions->True]

8. A bell-shaped plot for averaged (Subscript[M, coup])((Subscript[y, P]^( * )), (Subscript[y, Q]^( * )))

Frozen code

nn = 10;
k = 4;
Timing[Do[specs[coup] = Table[randomBoolNetOfCoupMULTIEDGE[nn, k, coup], {500}], {coup, 0, nn*k}]]

Timing[myData =
  Table[
   With[{mValsForCoup = Table[syncMutInfApprox[s], {s, specs[coup]}]},
    {coup, Mean[mValsForCoup], StandardDeviation[mValsForCoup]}],
   {coup, 0, nn*k}]]

ListPlot[Take[#,2]&/@myData,AxesLabel->{"coup","M((Subscript[Y, P]^( * )),(Subscript[Y, Q]^( * )))"},PlotLabel->Row[{500," (",nn,",",k,")","-nets per coup" }]]

Values of (Subscript[M, coup])((Subscript[y, P]^( * )), (Subscript[y, Q]^( * ))) averaged over 500 (10, 4)-nets per coup

(Subscript[M, coup=0])((Subscript[y, P]^( * )), (Subscript[y, Q]^( * ))) = 0.  Why?

(Subscript[M, coup=n*k])((Subscript[y, P]^( * )), (Subscript[y, Q]^( * ))) = 0.  Why?

Acknowledgements

This research and the associated interactive demonstrations were developed in the context of the FQXi-funded project ‘Mutual Information in weakly coupled process-algebraic mutual observers’ (grant FXQi-MGA-1702).
Thanks to Larissa Albantakis for exchanges on topics related with effective, mutual, integrated information and related concepts.

References

G. Tononi, An information integration theory of consciousness, BMC Neurosci. 5: 42, 2004.

D. Balduzzi, G. Tononi, Integrated Information in Discrete Dynamical Systems: Motivation and Theoretical Framework, PLoS Computational Biology, 4(6), e1000091, June 2008.

M. Oizumi, L. Albantakis, G. Tononi.  From the Phenomenology to the Mechanisms of Consciousness: Integrated Information Theory 3.0. PLoS Comput. Biol. 10(5), 2014, e1003588, doi:10.1371/journal.pcbi.1003588

L. Albantakis, G. Tononi, The Intrinsic Cause-Effect Power of Discrete Dynamical Systems -- From Elementary Cellular Automata to Adapting Animats, Entropy, 2015, 17, 5472-5502; doi:10.3390/e17085472.

S.A.Kauffman. Metabolic stability and epigenesis in randomly constructed genetic nets., J. Theoret. Biology, 22:437–467, 1969.

S. A. Kauffman, Homeostasis and Differentiation in Random Genetic Control Networks, Nature, 224, pp 177-178, 1969.

S. A. Kauffman, At home in the universe: the search for laws of self-organization and complexity. Oxford paperbacks. Oxford Univ. Press, 1995.

T. Bolognesi, LOTOS-Like Composition of Boolean Nets and Causal Set Construction, in: J-P Katoen, R. Langerak, A. Rensink (eds.), ModelEd, TestEd, TrustEd,  LNCS 10500 (Festschrift), pp. 27-47, Springer, 2017