Evolution enhances mutational robustness andsuppresses the emergence of a new phenotype– A study on gene regulatory network –Macoto Kikuchi (collaborator: Tadamune Kaneko)
Cybermedia Center, Osaka University

1

Contents

1. Introduction• robustness, gene regulatory network

2. Model and Methods• network model, multicanonical MC, evolutionary simulation

3. Results• ﬁtness landscape, bistability, mutational robustness

4. Discussions

2

Introduction

3

Central question
Mutational Robustness• Living systems do not easily lose their viability even when some of the genes
are mutated.
• This trait of mutational robustness has been developed through the long
history of evolution.
• Mutational robustness can be examined experimentally by several methods
such as gene knockout or artiﬁcial addition of gene regulation.
• e.g. Isaran et al. 2008.

• But the answer to the question, “Does mutational robustness evolves?”
remains unclear because experimental investigations are diﬃcult.
• Thus numerical methods play an indispensable role.
4New methods are required
We investegated the evolution of mutational robustness for gene regulatory
network (GRN) model.
Problem• Evolutionary simulations only provide evolutionary-path-dependent traits.
• To study the universality and the characteristics of evolution, we need areference system.
• The reference system that we consider appropriate is a set of randomly
sampled GRNs.
• A simple random sampling is useless because highly ﬁt GRNs are rare.
5

McMC to study evolution
Solution• Rare event sampling using the multicanonical Monte Carlo (McMC) method• Random sampling over wide range of ﬁtness

The preceding works applying McMC to GRN• N. Saito and M. Kikuchi: New J Phys 15 (2013) 053037• S. Nagata and M. Kikuchi: PLoS Comput Biol 16 (2020) e1007969
This talk• T. Kaneko and M. Kikuchi: arXiv:2012.03030
6

Gene Regulatory Network 1
Gene expression• Information coded in agene is transcribed tomRNA.
• Ribosome assembles
amino acids according tomRNA sequence.
• Protein is produced.

7

Gene Regulatory Network 2
Gene regulation• A category of proteins
calls transcription facotr
regulate expression ofother genes working either
as activator or repressor.

8Gene Regulatory Network 3
Gene regulatory network• Genes regulate mutuallyand form a complex
network to controle thecell states.
• Adaptation toenvironmental change• Maturation,
diﬀerentiation

9Model and Methods

10Model 1
Connectionist model• We ignore details of gene expression andtake into account only regulatoryinteractions.
• Node: Gene
• Edge: Regulation

Details• GRNs are represented by directed random graphs• N = 32 nodes、K = 80 edges• One input gene and one output gene
• All the possible networks are allowed11Model 2
Discrete-time dynamics(xi (t + 1) = RI δi,0 +∑)Jij xj (t)jR(x) =

1
1 + e −β(x−µ)• 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])• µ = 0, β = 2 (Hopﬁeld-Tank)

12

Model 3Fitness• We require that the expression of the output node for I = 0 and 1 diﬀer aslarge as possible.
• Fitness is deﬁned as follows:
f = |x̄out (0) − x̄out (1)|• x̄out (I ) : the ﬁxed-point value (or the long-time average if no ﬁxed point isreached) of the output node for input I .

• f ∈ [0, 1] by deﬁnition
13

Methods 1

Two methods
1. Random sampling by multicatnonical Monte Carlo method
2. Evolutionary simulation

14

Method 2
Multicanonical Monte Carlo• Originally developed in statistical physics to sample energies evenly.
• Berg and Neuhaus 1991,1992• Later, it was realized that McMC can used also for non-physical systems.
• Review: Y. Iba, N. Saito, and A. Kitajima: Ann Inst Stat Math 66 (2014)
611.

Application to GRN• It enables us to sample GRNs in wide range of ﬁtness randomly (in principle).
• Fitness are divided into 100 bins.
• The weight for each bin is determined by Wang-Landau method.
• 5M samples

15

Method 3
Example of McMC
Histgram

Probability
050000−2−4logΩ(f)

40000

30000

20000

−6−8−10−12−14

10000

−16
0
0.0

0.2

0.4

0.6

f

0.8

1.0

0.0

0.2

0.4

0.6

0.8

1.0

f16

Method 4

Evolutionary simulation• Mutation
1. Delete a randomly selected edge
2. Add a new egde to a randomly selected node pair

• Initial population: 1000 random GRNs.
• 500 are selected at each generation based on ﬁtness.
• The lineage of the highest ﬁtness at 150th generation is sampled.
• 100,000 independent runs: 100,000 lineages

17

Results

18

Fitness landscape and evolutionary speed 10

Ω(f)

4060
random sampling
evolution

0.0

0.2

0.4

0.6

0.8

1.0

100

0

10−1

20

10−2

40

10−3

80

10−4

100

10−5f

Generation

Generation

20

Ω(f)

10−1
10−3
10−5
10−7
10−9
10−11
10−13
10−1560
random sampling
evolution

80
1000.00.2

0.4

0.6

0.8

1.0f

Probability of ﬁtness and ﬁtness vs.
generation for β = 2for β = 4 (only this ﬁgure)
19

Fitness landscape and evolutionary speed 2Observations
1. Most GRNs have low ﬁtness (f < 0.1).
2. At intermediate ﬁtness, GRNs become exponentially rare as ﬁtness increases.
3. For very high ﬁtness (f > f ∗ ), GRNs become faster-than-exponentially rare.
4. Evolution slows down substantially at f ≃ f ∗ .
Genotipic entropy matters• The evolutionary speed is determined mainly by the number of available
GRNs (or genotypic entropy).
20

Emergence of bistability 1

1.0

1.0
up

random samplingevolution

down

0.8

0.8x(I)P2 ( f )

0.6

0.4

0.6
0.4
0.2

0.2

0.0
0.0
0.0

0.2

0.4

0.6

I

0.8

1.0

0.2

0.3

0.4

0.5

0.6

f

0.7

0.8

0.9

1.0

An example of bistability for f = 0.9 The fraction of the bistable GRNs
21

Emergence of bistability 2
Universality of evolution• Bistable GRNs emerge, and their fraction grows as ﬁtness increases. AllGRNs become bistable as f → 1.
• Bistability is a new phenotype. This new phenotype that was not explicitly
required by the ﬁtness appears inevitably.

Particularity of evolution• Evolution delays the appearance and growth of the fraction of bistable GRNs
Evolution is conservative!
22

Mutational robustness 1Robustness measure• We consider the single-edge deletion as mutation.
• Delete a single edge from GRN and compute ﬁtness f ′ after the mutation.
• Robustness measure for a given GRN is the average of f ′ for all the possible
single-edge deletion.
1 ∑ ′r≡fK i i• fi ′ is f ′ for GRN that ith edge is deleted.

23

Mutational robustness 2• Since ⟨r ⟩ increases with f , it is not suitable to compare the robustness
measure of diﬀerent ﬁtness f .
• We focus on diﬀerence in ⟨r ⟩ of GRNs having the same f obtained byrandom sampling and evolution.

24

Mutational robustness 3
Robustness measure
1.0

r⟩ndom⟨s⟩mpling
evolutionTwo-step evolution

0.8

1. Early: ⟨r ⟩ coincides withrandom sampling.

⟨r⟩

0.6
0.4

2. Later: ⟨r ⟩ exceeds random
sampling.

0.2
0.0
0.0

0.2

0.4

0.6

0.81.0f25

Mutational robustness 4
Probability distribution of r
0.200

random sampling
evolution

0.35
0.30

P( r )

0.15

P( r )

0.125

0.20

0.100
0.075

0.10

0.050

0.05

0.025

0.0

0.2

0.4

r

0.6

0.8

1.0

0.06
0.04
0.02

0.000

0.00

random sampling
evolution

0.08

0.150

0.25

P( r )

0.10

random sampling
evolution

0.175

0.00
0.0

0.2

0.4

r

0.6

f = 0.8
f = 0.5
Evolution enhances mutational robustness!

0.8

1.0

0.0

0.2

0.4

r

0.6

0.8

1.0

f = 0.99

26

Mutational robustness 5
Lethal edges
0.5

P(f')

0.4

• Edges are divided in twoclasses: neutral and lethal

0.3

0.2

• Intermediate edges arescarece.

0.1

0.0

0.0

0.2

0.4

f'

0.6

0.8

Distribution of f ′ for f = 0.99

1.027

Mutational robustness 6
Distribution of lethal edges
random sampling
evolution

0.08

P(nL)

0.06
0.04
0.02
0.00
010

20

30

40

nL

50

60

7080

• Lethal edges are
signiﬁcantly scarce inevolutionally obtainedGRNs.
• The origin of thediﬀerence in themutational robustness
is the diﬀerence in thenumber of lethal edges.
28

f = 0.99

Mutational robustness 7Eﬀective network size
0.200

random network
random sampling
evolution

0.175
0.150

• We counted the number ofnodes having at least one
path to the output node.

P(nN)

0.125
0.100
0.075

• The eﬀective network
size does not aﬀect the
mutational robustness.

0.050
0.025
0.000
051015nN

20253029

f = 0.99

Network motifs 1
Deﬁnition of motifs• The connection patterns appear in the given network at frequencies muchhigher than those found in random networks.

• We counted the following patterns for f = 0.99
1. auto-regulation
2. mutual regulation
3. feedforward loop
4. feedback loop
5. mutual activation/repression plus auto-activation
30

Network motifs 2

0.5

random_network
random_sampling
evolution

0.30

0.4

P(nm + )P(nS + )

0.25
0.20
0.15
0.10

0.3
0.2

0.00
246

nS +

8

Auto-activation

10

12

0.4
0.3
0.2

0.0

0.0
0

random network
random sampling
evolution

0.1

0.1

0.05

0.5

random network
random sampling
evolutionP(nm − )0.350246nm +810

Mutual activation

120246

nm −

81012

Mutual repressionn

31

Network motifs 3

0.150

0.30

random network
random sampling
evolution

0.100
0.075
0.050

0.20
0.15
0.10

random network
random sampling
evolution

0.8

P(nMA − )

0.25

0.125

P(nFBL + )P(nFFL + )1.0

random network
random sampling
evolution

0.175

0.6
0.4
0.2

0.025

0.05

0.000
0510

nFFL +

1520

Coherent feed-forward
loop

0.0
0.00

00246

nFBL +

810

Positive feedback loop1212

nMA −345mutual repression plus
auto-activation

32

Network motifs 4• Six patterns may be called as motif (but not so signiﬁcant)• +FFL is ubiquitous motif in GRNs
• Other patterns are known as motifs in multistable GRNs.

• Motifs are related to the function but not to the mutational robustness.
• Global structures are considered to be relevant to the mutational robustness.

33

Path distribution 1
The number of paths connecting the input
and output nodes
0.006

random sampling
evolution

0.005

P(npath)0.004• npath are signiﬁcantly
diﬀerent for a small
number of paths

0.003
0.0020.001• npath = 1 in particular.

0.000−0.001050

100

150200npath

250300

350

400

34

Path distribution 2
The number of paths vs. the number oflethal edges• The lethal edges are
signiﬁcantly scarce inevolution, irrespective of
the number of paths.
• The diﬀerence npath
does not fully explain
the mechanism forenhancement of themutational robustness.
35

Discussions

36

Mechanism for enhancement of mutational robustness
Mutation consists of two successive processes:
1. Deletion of a randomly selected edge
2. Addition of a new edge
• When a lethal edge is deleted in the ﬁrst process, the ﬁtness of the GRN
drops down to almost zero, and the possibility is very low that the ﬁtness
recovers by the second process.
Second-order selection (Wagner)• GRNs with many lethal edges are more likely to be eliminated than thosewith fewer lethal edges in evolution.
37

Evolutionary speed
Entropy eﬀect• Speed of evolution is in large part determined by the number of available
GRNs.
• Eﬀect of genotypic entropy.

Question• Sato et al.(2003) showed that evolution slowed down as ﬁtness becomeshigh.
• Consistent with the entropic eﬀect• They found that evolution speed is in association with phenotypicdivergence.
• Relationship between the phenotypic divergence and the genotypic entropy?
38

Summary39

New method• We compared the evolution of GRNs with a randomly sampled set generated
by the multicanonical Monte Carlo.
Results
1. Mutational robustness does evolve.
• GRNs with fewer lethal edges are selected by evolution.
• Motifs are not relevant to mutational robustness.
• The number of paths has some relevance to mutational robustness.

2. Evolution delays the emergence of a new phenotype.
• Evolution is conservative.

3. Evolutionary speed is determined mainly by genotypic entropy.
40