Function and Robustness of GeneRegulatory Network: Toward theLandscape Picture of Evolution
Macoto Kikuchi (Collaborators: S. Nagata andT. Kaneko)
大阪大学サイバーメディアセンター

2019/3/7

Contents

1
2
3
45

Motivation
Gene Regulatory Network(GRN)
ModelMethod (Multicanonical Monte Carlo)
Results: Function and Robustness

Motivation

Characteristic properties of ”evolved
thing” are function and robustness
A.Wagner: ”Robustness and Evolvability in Living
Systems” (2005)
Intuitively, highly optimized system may befragile.
Evolution is not simply an optimization process?

Robustness
Robustness against perturbation
Stability in development and diﬀerentiation:Canalization (Waddington)
Epigenetic landscapeProtein folding: Anfinsen’s dogma, Funnel picture(Go, Wolyness)
Energy landscape

Robustness against mutation
Function is not lost by mutation
Homologous protein

Prospect
Landscape picture of evolution
Consider evolution landscape, including phenotypes
not visited in the course of evolution
Evolutional pathway on the landscape

We consider a toy model of the gene regulatory
networkAs the evolved system should be rare, we use
the rare event sampling method

Gene Regulatory Networks(GRN)Gene expression

Gene regulation

A complexnetwork in which
the genes mutually
regulate by thetranscription
factors (TF)
Abstract model of GRN

TFs themselves
are proteinsmade by thegene expressions

Question
Character of the fitness landscape
Relation between the cooperative response tooutside and the robustness
Mutational robustness
Robustness against external/internal fluctuation(number fluctuation of TF or other molecules)

Can we see any universal properties, if we classify
the randomly generated GRNs by fitness
Properties that do not depend on theevolutional pathway

ModelSimple toy model of GRN having one input geneand one output genecf. M. Inoue and K. Kaneko PLOS Compt.
Bio. 9(2013)e1003001
Directed random graph: N nodes、K edges
Node: Gene
Edge: Regulatory relation
Self regulation and mutual regulation are
excluded
The input node is randomly selected from thenodes having paths to all the other nodes
The output node is selected from the nodes

IO
GRN having one input node and one output nodeand having no self and mutual regulations

Discrete-time dynamics (Neural-network like)Xi (t + 1; I ) = R (I δj,1 + Σj Jij Xj (t; I ))R(x) =

tanh x + 1
2

cf. A. Wagner: Evolution 50 (1996) 1008
Xi : Expression of ith gene ( [0, 1])
Jij : Regulation of ith gene by jth gene(0, ±1)+1: activation, −1: repression

I : Input from exterior world ([0, 1])
R(x) : Soft response function

1.0

0.8

R(x)

0.6

0.4

0.2

0.0
−2.0

−1.5

−1.0

−0.5

0.0

x

0.5

1.0

1.5

Response function

2.0

Spontaneour
expression is 0.5(comparativelylarge)
We want toassemble a circuitthat can respondsensitively to Onand Oﬀ ofexternal signal
cf. M. Inoue andK. Kaneko:EPL124 (2018)
38002

Required function
Sensitive response to On-Oﬀ changeof external signal
Since the response damps out for sequential
circuit, Feed-Forward type regulation isindispensable
Both activation and repression are requiredFitnessX̄i (I ): Temporal average of the response of ithgene (in the steady state)
Sensitivity of ith gene: Deﬀerence of theresponse to I = 0 and 1
Si = |X̄i (1) − X̄i (0)|The node having the largest sensitivity isdefined as the output node
Xout : X of the output node (response of thenetwork)
Fitnessf ≡ S : Sensitivity of the output nodeMethod (Multicanonical MC)

Sampling methodthat gives a flat
distribution ofenergy

Ideal energy histgram
obtained by themulticanonical MC

Enable us tosamplelow-energy rarestatesEnable us tocalculate thedensity of states

Detailed balance
wij P(Ej ) = wji P(Ei )For ordinary Metropolis MC, P(E ) ∝ e −βE

We can use any P(E )P(E ) ∝ e −f (E )and we requiree −f (E ) ∼1Ω(E ):Weight f (E ) is determined through learningprocess

Using the obtained energy histgram, we canestimate DOSΩ(E ) ∝ H(E )e f (E )Divide E into bins
Piecewise linear approx. for f (E ): Multicanonical
B.A. Berg and T. Neuhaus: PRL 68 (1992) 9Constant f (E ) in each bin: Entropic sampling
J. Lee: PRL 71 (1993) 211

Wang-Landau method for the learning processused only for the entropic sampling
F. Wang and D.P. Landau: PRL 86 (2001) 2050

2D Ising Model

6000

800
5000600400

4000

EH(E)

200300002000−200−4001000−6000−8000

200000400000MCS

6000008000001000000

Time series of energy−800−600−400−2000E

200

400600800

Obtained energy
distribution

250logΩ(E)200
150100500−800−600−400−2000E

200

400600800

DOS estimated bymulticanonical MC

Number of theground state is
2.07 (cf. truevalue is 2)

Application to non-energetic system
Eigenvalue distribution of random matrix
N. Saito, Y. Iba and K. Hukushima: PRE 82(2010) 031142

Search for periodic orbits in a chaotic system
A. Kitajima and Y. Iba: Compt. Phys. Comm.
182 (2011) 251Stability of a coupled chaotic map
N. Saito and M. Kikuchi: New J. Phys. 15 (2013)053037

Enumeration of magic squares
A. Kitajima and M. Kikuchi: PLOS One 10 (2015)
e0125062

The first paper that discuss the evolutionary
landscape using multicanonical MC”Robustness leads close to the edge of chaos incoupled map networks: toward the understanding ofbiological networks”
N. Saito and M. Kikuchi: New J. Phys. 15 (2013)053037
Evolution and robustness of a coupled chaotic map(an abstract model for GRN)

Application to GRN
Sampling that gives the flat distribution offitness
Divide the fitness (0 ∼ 1) into 100 bins

In principle, we can randomly sample GRNs
with several diﬀerent values of fitness
Actually there are correlations between samples

Microcanonical ensemble within each binN = 16 ∼ 32
Average number of edges connected to eachnodeC ≡ 2N/K = 5, 6
We show results for C = 5 mainly

Results 1

Fitness Landscape(a)C = 5(b)C = 6
N=32 K=80
N=28 K=70
N=24 K=60
N=20 K=50
N=16 K=40

10−2
10−510−5P( f )10−8P( f )10−810−1110−1110−1410−1410−1710−1710−20

N=32 K=96
N=28 K=84
N=24 K=72
N=20 K=60
N=16 K=48

10−2

0.00.20.4

f0.6

0.8

1.0

10−200.0

0.20.4

f0.6

0.8

1.0

Appearance probability Appearance probabilityof fitness (C = 5)of fitness (C = 6)

There is a threshold of rarenessmore than 95% in f < 0.2

GRNs having f larger than the threshold areexponentially raref > 0.9 are more than exponentially raref > 0.99: The fittest ensemble

GRNs with high fitness are rareResponse in steady states(a)f = [0.7, 0.71]
1.0

Steady-state
response when theinitial condition is
Si = 0.5 for all ixout I

̄ ()

0.8

0.6

0.4

0.2

0.0
0.0

0.2

0.4

0.6

0.8

I

f = [0.7, 0.71] (20samples)

1.0

Smoothresponse to the
input IA single fixedpoint

(b)f = [0.99, 1]
1.0

xout Ī ()0.8

Step-like responseto the input I

0.6

0.4

0.2

0.0
0.0

0.2

0.4

0.6

0.8

1.0I

The fittest ensemble(20 samples)

Response byswitching two
fixed points
Ultrasensitivity

Responses for f ∼ 0.7(b)1.0(a)
1.0

forwardbackward(1000; I)

0.8

0.6

xout

xout

(1000; I)

0.8

0.4

forwardbackward

0.6

0.4

0.2
0.2
0.0
0.0

0.0
0.0

0.2

0.4

0.6

0.8

1.0

0.2

0.4

0.6

0.8

I

I

Sweeping I (nobifurcation case)

Sweeping I(saddle-nodebifurcation case)

1.0

Responses of the fittest ensemble(c)1.0

forwardbackward

xout(1000; I)

0.8

0.6

0.4

0.2

0.0
0.0

0.2

0.4

0.6

0.8

1.0

I

Sweeping I (saddle-node bifurcation)

Appearance probability of two fixed points
Monotone increaseagainst the fitness
1.0Correspondence
between thefunction and thenumber of thefixed points

N=162024
28
32

0.8

P2

0.6
0.4
0.2
0.0
0.5

0.6

0.7

f0.8

0.9

1.0

99% of GRNs in
the fittest
ensemble have two
fixed points

As the fitness increases, the big jump that the
number of the fixed points changes takes place atsomewhere in the course of evolution, irrespective of
the evolutionary pathway
Universality of evolution
The fitness restricts the phenotypedynamical response
1.0

Xout (t; I)

0.8

61% of the GRNs
can respondproperly.

0.6

0.4

0.2

0.0
0200

400600

800

1000

1200

1400

t

Dynamical response tosudden changes of theinput

Whether or notthe bistable
range include 0or 1

Robustness against the input noise
Numberfluctuation of theinput molecule
I

1.2

xout

1.0

0.8

0.6

0.4

0.2

0.0

−0.2

0200

400

600

800

t

1000

1200

1400

Uniform noise of[−0.3, 0.3]

GRNs that canrespond to thesudden change ofinput are robust
against the input
noise
The eﬀect ofthe fixed-pointswitchingRobustness against the internal noise
Numberfluctuation of TF

(a)
1.0

without noisewith noise

xout t I(;)

0.8

0.6

0.4

0.2

0.0
0200

400

600

800

1000

1200

1400t

Dynamical responsewhen the internal noise
is applied

Uniform noisesof [−0.2, 0.2]
are applied to all
the input toeach gene

GRNs that can
respond to thesudden change ofinput are robust
against theinternal noise

Noise-induced ultra sensitivity(b)0.80.8(;)

1.00.6

xout t I

xout t I

(;)

(b)
1.0

0.4

0.20.60.40.2
without noisewith noise

0.00200400

without noisewith noise

0.0
600

800t

1000

1200

14000200400600800t

1000

1200

1400

Fixed points and the robustness

60% of the GRNs having two fixed-points can
respond properly to the sudden change of input
They are robust against both input and internal
noisesSome of the GRNs exhibit the noise-induced
ultrasensitivity
70% in total can respond properly to the inputMutational robustness

Mutation of the single-edge deletion
A moderate mutation (e.g. slight change of TF)
We try all the possible mutations
Input/Output nodes are unchanged upon mutation

Distribution of the fitness f ′ after themutation
(b)

(a)
0.9
0.8
0.7
0.6
0.5

0.40
0.35

0.35
0.30
0.25

0.25′

P( f )′P( f )

0.30

0.20

0.20
0.15

0.15

0.10

0.10
0.05

0.05
0.00

0.99

0.40

0.0

0.2

0.4

f

0.6

0.8′

Several diﬀerent f

1.0

0.00

0.0

0.2

0.4

f

0.6

0.8

1.0′

The fittest ensemble

Majority of the edges are neutral against
mutation
For the fittest ensemble, most of the edges areeither neutral or lethalIntermediate edges are scarse

(a) C = 5
N=1620
24
28
32

0.10
0.08

P(nL)

0.06
0.04
0.02
0.00
02040

nL60

80

100

Distribution of lethal
edges

Typical number ofthe lethal edges is6
Independent ofsize
Larger GRNs are
relatively robust
Some GRNs
have no lethaledge

IOAn example of GRN without a lethal edge

Results 2

Compare the results of evolutionary simulations
and the random sampling
Slightly diﬀerent model just for simplicity
Allow both self and mutual regulation. Fixed
input/output nodes

Population: 1000
Keep 500 samples from the highest fitness. Applymutation to 500 copies
Perform 10000 runs and follow the evolutionary
path of the fittest sample

Generation and fitness Fitness landscape and(fitting by tanhthe speed of evolution

Speed of evolution is determined by entropyRobustness index
For all the possible deletion of edges,
1 ∑ ′fKedge

is defined as the robustness of each GRN

Robustness distribution
Avarage robustness andand the evolutionary
evolutionary pathway
pathwayDistribution of therobustness for thefittest ensemble
Distribution of therobustness forsamples just afterf > 0.99 is
attained for theevolutionary
simulations.

Evolutionary process is divided into two stages
1
Entropic stage
2
Robustness-aquiring stage

Pathways for diﬀerent
number of copies

Summary

GRNs of high fitness have the following features:
1
Ultrasensitivity (two stable fixed points)
A big jump irrespective of the evolutionary pathway
2

Three robustnesses
mutation
input noiseinternal noise

Evolution enhances robustness