  The aim of this paper is two-fold. First, we propose a new computational
method to investigate the particularities of evolution. Second, we apply this
method to a model of gene regulatory networks (GRNs) and explore the evolution
of mutational robustness and bistability. Living systems have developed their
functions through evolutionary processes. To understand the particularities of
this process theoretically, evolutionary simulation (ES) alone is insufficient
because the outcomes of ES depend on evolutionary pathways. We need a reference
system for comparison. An appropriate reference system for this purpose is an
ensemble of the randomly sampled genotypes. However, generating high-fitness
genotypes by simple random sampling is difficult because such genotypes are
rare. In this study, we used the multicanonical Monte Carlo method developed in
statistical physics to construct a reference ensemble of GRNs and compared it
with the outcomes of ES. We obtained the following results. First, mutational
robustness was significantly higher in ES than in the reference ensemble at the
same fitness level. Second, the emergence of a new phenotype, bistability, was
delayed in evolution. Third, the bistable group of GRNs contains many
mutationally fragile GRNs compared with those in the non-bistable group. This
suggests that the delayed emergence of bistability is a consequence of the
mutation-selection mechanism.
