In the cytoplasm of a biological cell there are a myriad of di erent proteins, lipids and enzymes, each of them performing di erent tasks. Therefore, the spatial organisation of these chemical species is crucial for the functioning of a cell. Order can emerge in many ways but, in living beings, it often appears as a consequence of free-energy expenditure. Hence, in this thesis, we study how order and organisation can emerge from the interplay between physical interactions and free-energy consumption, both analysing particular examples that arise in cell biology and constructing a more abstract framework for reaction-di usion systems in the presence of interactions.

In the rst part of this thesis, we explore, from a theoretical perspective, two processes where the intracellular medium is patterned and organised, whose common feature is the fact that they both stem from the collective behaviour of a large number of molecules. First, we develop a model for protein aggregation which studies the e ect of intracellular obstacles on the coagulation kinetics. Our predictions are then successfully compared with experimental data.

The other example refers to the compaction and location of the bacterial chromosome. We suggest that the chromosome segregates from the rest of the cytoplasm, because of steric interactions between DNA and the intracellular crowders, by means of a mechanism reminiscent of liquid-liquid phase separation. Also, our study indicates that spatial localisation within the cell is dictated by non-equilibrium transcription of mR-NAs (which are part of the crowding e ect). Our model successfully reproduces much of the phenomenology observed in Escherichia coli cells, for example, the positioning of the chromosome at di erent stages of the cell cycle.

Building on these examples, in the second part of the thesis, we construct a thermodynamically consistent framework to mathematically describe chemical reaction networks in non-ideal solutions. This framework allows us to generalise the results from the classical theory of ideal networks, and aids in elucidating the connection between non-equilibrium chemical reactions and phase separation for a large class of networks, known as complex-balanced networks. Given that complex balance can be fully determined from the topology of the network, we analyse how this topological property of the network can constrain the dynamics of the solution, and what behaviour we can expect when complex balance is broken.

We conclude by discussing how the two approaches used here are related and the contribution each of them can have for the advancement of biophysics.

Résumé

La cellule est l'élément de base composant tous les êtres vivants sur Terre. À l'intérieur de chaque cellule, dans le cytoplasme, il y a une myriade de constituants tels que des protéines, des lipides et des enzymes. Chacune de ces molécules joue un rôle bien spéci que au sein de la cellule, ce qui requiert un contrôle précis de leur position spatiale. L'organisation intracellulaire est donc essentielle.

Pour parvenir à cette organisation spatiale détaillée, les êtres vivants ont recours à stratégies diverses. D'un côté, la mise à pro t des interactions entre les constituants de la cellule peut créer de l'ordre sans nécessiter de dépense énergétique, mais cela ne peut créer qu'une organisation statique. D'un autre côté, il existe d'autres stratégies qui sont plus dynamiques mais qui nécessitent une dépense d'énergie (techniquement, une dépense d'énergie libre) qui peut, parfois, être très importante.

L'une des caractéristiques des organismes vivants est la dépense continue d'énergie. Par conséquence, pour comprendre comment l'auto-organisation des organismes vivants émerge, il est nécessaire de comprendre comment cette dépense énergétique modi e les lois physiques que nous connaissons pour les systèmes à l'équilibre thermodynamique.

Il y a donc deux grands axes de recherche dans cette thèse, chacun d'entre eux étant représenté dans une des parties de cette thèse. Dans la première partie, nous proposons une explication théorique pour deux processus où le milieu intracellulaire est organisé spatialement sous forme de motifs. Ces processus découlent tous deux du comportement collectif des molécules. Dans la deuxième partie, nous construisons un cadre thermodynamique pour décrire les réseaux de réactions chimiques dans des solutions non idéales, c'est-à-dire des solutions avec des interactions entre solutés. Cette approche nous permet d'obtenir des résultats très généraux pour les modèles de solutions non idéales avec des réactions chimiques hors équilibre thermodynamique. Cette classe de modèles peut décrire une grande variété de phénomènes dans les cellules biologiques.

Première partie

Cette partie est consacrée à l'analyse théorique des deux processus au cours desquels le milieu intracellulaire s'organise et des motifs spatiaux peuvent apparaître. Dans un premier temps, nous développons un modèle qui a pour but d'étudier l'impact de l'encombrement stérique intracellulaire sur la cinétique d'agrégation des protéines (Chapitre 1). Dans un deuxième temps (Chapitre 2), nous étudions la compaction et la localisation spatiale du chromosome de la bactérie Escherichia coli (E. coli).

iii

Encombrement stérique et agrégation de protéines

Dans l'environnement intracellulaire, le cytosquelette et d'autres grands obstacles entravent considérablement la di usion des particules. Cependant, la manière dont ces obstacles a ectent les particules dépend de la taille de chacune d'entre elles. Certains obstacles, comme le cytosquelette, n'a ectent que les particules de plus de 35-50 nm [START_REF] Luby-Phelps | Hindered di usion of inert tracer particles in the cytoplasm of mouse 3T3 cells[END_REF], Etoc et al., 2018], un seuil parfois appelé "taille des pores du cytoplasme". Pour ce type de particule (taille d'environ 50nm ou plus), le cytosquelette ou d'autres obstacles de grande taille ralentissent fortement la di usion des particules et peuvent a ecter de nombreux processus ; c'est le cas de l'agrégation des protéines.

Nous étudions, d'un point de vue théorique mais complété par des expériences, comment ce type d'encombrement stérique modi e la cinétique d'agrégation des protéines à l'intérieur de la cellule. Nous développons un modèle mathématique pour rendre compte de ces deux coe cients de di usion très di érents et, par comparaison avec les données expérimentales, nous trouvons comment l'encombrement stérique conditionne le nombre et la taille des agrégats des protéines. En n, nous calculons approximativement la taille critique des agrégats au-delà de laquelle ils sont piégés dans le cytoplasme cellulaire, et nous constatons qu'elle est proche de la taille du pore du cytoplasme, ce qui est en accord avec notre modèle.

Le nucléoïde chez E. coli

L'organisation intracellulaire des bactéries est un cas d'étude très intéressant, car aucune membrane ne sépare les di érentes parties de la cellule. Cependant, les bactéries présentent un très haut degré d'organisation intracellulaire. Le cas du chromosome de la bactérie E. coli, qui est situé dans une région au centre de la cellule que nous appelons le nucléoïde, est particulièrement spectaculaire. En outre, le nucléoïde contrôle et organise spatialement de nombreux processus cellulaires. Par exemple, il découple spatialement la transcription et la traduction en protéines car les ribosomes, lorsqu'ils sont xés à l'ARN messager (ARNm), sont exclus du nucléoïde à cause des interactions stériques [START_REF] Sanamrad | Single-particle tracking reveals that free ribosomal subunits are not excluded from the Escherichia coli nucleoid[END_REF]. Finalement, la position du nucléoïde change beaucoup au cours du cycle cellulaire : Quand il n'y a qu'un seul chromosome, il est situé très précisément à la moitié de la cellule, mais juste avant la division cellulaire, les deux chromosomes se positionnent à un quart et trois quarts de la cellule (le long de l'axe principal de la cellule) [Wu et al., 2019].

Pour tenter d'expliquer ces phénomènes, nous avons développé un modèle qui prend en compte les interactions stériques entre le chromosome et les «crowders» de grande taille comme les ribosomes ou les polysomes (ARNm avec un ou plusieurs ribosomes attachés). Même s'il y a d'autres e ets a ectant potentiellement la taille et l'emplacement du nucléoïde, notre hypothèse de base est que, pour un grand nombre des cas expérimentaux, l'e et dominant est l'encombrement stérique. Notre modèle montre que iv l'encombrement stérique (couplé à des e ets intrinsèquement hors équilibre thermodynamique, comme la transcription et la traduction en protéines) peut expliquer l'apparition du nucléoïde (par un mécanisme similaire à celui de la démixtion de deux liquides). Le modèle peut également expliquer pourquoi la taille de nucléoïde augmente lorsque la transcription d'ARNm est arrêtée et diminue lorsque la traduction en protéines est arrêtée. D'autre part, selon notre modèle, le positionnement central du chromosome, et à 1/4 et 3/4 pour cellules plus grandes, est une conséquence directe de la transcription des ARNm: les ARNm participent à l'encombrement stérique et donc altèrent la pression osmotique dans les di érentes zones de la cellule, créant des forces qui pussent le nucléoïde au centre (ou à 1/4 et 3/4, en fonction de la taille du cellule).

En n, comme tous les autres modèles, le nôtre présente également des limites. L'une de ces limites est l'hypothèse d'homogénéité dans la direction radiale de la cellule, ce qui implique que le système soit e ectivement unidimensionnel et facilite les calculs. Néanmoins, dans certains cas, l'homogénéité radiale est rompue (notamment si le nucléoïde se contracte radialement), ce qui peut avoir des conséquences importantes sur le positionnement du nucléoïde. Par conséquent, nous prédisons, comme conséquence de notre modèle, que si une perturbation contracte le nucléoïde radialement, elle déstabilise également son positionnement.

Deuxième partie

Dans l'analyse précédente, nous avons considéré quelques exemples de phénomènes au sein des cellules qui ont pour point commun de provenir des aspects non idéaux des solutions, notamment l'encombrement stérique et la démixtion des liquides. D'autre part, nous avons vu que des phénomènes hors équilibre thermodynamique (particulièrement des réactions chimiques) contribuent à l'organisation cellulaire. Donc, dans la deuxième partie de cette thèse, nous construisons un cadre thermodynamique général pour décrire les réseaux de réactions chimiques au sein des solutions non idéales.

Cette approche, à la fois générale et thermodynamiquement cohérente, nous permet de généraliser certains résultats de la théorie des réseaux de réactions chimiques idéales tels que les fonctions de Lyapunov (fonctions qui sont minimisées par la dynamique du système) pour les réseaux de type «complex balanced» [START_REF] Horn | General mass action kinetics[END_REF]. Ce dernier type de réseau o re un cadre moins restrictif que les réseaux de type «detailed balanced» (qui sont des systèmes à l'équilibre thermodynamique) et, selon la théorie classique de réseaux de réactions chimiques, les réseaux de type «complex balanced» peuvent être caractérisés par leur topologie [Feinberg, 1972]. D'autre part, ce cadre nous permet de formaliser des liens entre les réactions chimiques hors équilibre et les séparations de phase (ou démixtions) pour des réseaux de ce type, pour lesquelles nous pouvons obtenir des diagrammes de phase complets en minimisant leur fonction de Lyapunov. Il est important de noter que cette classe de réseaux peut être identi ée uniquement à partir de la v topologie des réseaux, et que l'existence d'une fonction de Lyapunov est indépendante de la cinétique particulière ou de l'énergie libre de chaque système.

Finalement, nous montrons que cette topologie peut contraindre la dynamique des solutions et leur stabilité. En fait, nous montrons que, pour les réseaux de type «complex balanced» avec des interactions locales, l'apparition de «non-equilibrium patterning» (motifs dus à des réactions chimiques hors équilibre) est interdite. À l'inverse, notre théorie prédit que lorsque le réseau n'est plus de type «complex balanced», des dynamiques plus exotiques peuvent apparaître, comme celles de type «Ostwald ripening arrest».

Conclusion

En résumé, cette thèse a présenté deux approches di érentes de la même problématique, à savoir l'organisation intracellulaire en présence d'e ets hors équilibre thermodynamique. Dans la première partie, nous avons étudié quelques exemples d'un point de vue plutôt e ectif, ce qui nous a permis d'analyser de manière simple quelques processus cellulaires et de faire des prédictions qui peuvent être validées expérimentalement. Dans la deuxième partie, en revanche, nous avons utilisé une approche plus générale et abstraite, pour étudier le comportement d'une très grande classe de modèles qui pourraient être utiles pour la modélisation de l'auto-organisation cellulaire. Cependant, cette abstraction a rendu plus di cile la validation de la théorie et nous n'avons pas réussi à décrire certains phénomènes très complexes comme celui du «non-equilibrium patterning». Néanmoins, pour faire progresser la biophysique théorique, une intégration de ces deux approches sera nécessaire, ainsi qu'une synergie réussie entre théorie et expériences.

vi First and foremost, I would like to thank Michele for giving me this opportunity, for your time, and your patience. I would also like to thank Jean-François and Ned for the discussions we have had, always insightful. Also, collaborating with Mathieu and Aléria was one of the great pleasures of the early days of my PhD, for which I am grateful and from whom I learnt a lot. Finally, I would like to thank the rapporteurs, John Palmeri and Ulrich Gerland, who have agreed to evaluate my work, as well as the rest of the members of the jury, David Lacoste and Andela Šarić. Un grand merci aussi à Romain et Thomas, car si le texte français de cette thèse n'est pas trop mauvais c'est surtout grâce à vous.

On a more personal note, I would like to thank friends and colleagues from the unit (in fact, the UMR168 as a whole), for a great atmosphere for research but also for your warm welcome and making me feel as one more of the group from the very rst day. It has been a pleasure to do research here and it has mostly been so because of you. A special thanks to the people in the theory group (past and present) and the people that have passed through the o ce, for the always stimulating discussion and company, about science and beyond science.

Y, nalmente, un inmenso gracias a mi familia, por todo su apoyo durante estos años de duro trabajo.

P

Science is an essentially anarchic enterprise: theoretical anarchism is more humanitarian and more likely to encourage progress than its law-and-order alternatives. [...] And is it not clear that successful participation in a process of this kind is possible only for a ruthless opportunist who is not tied to any particular philosophy and who adopts whatever procedure seems to t the occasion?

This apparent paradox is resolved by noting that, while the entropy of the entire universe must always increase, it may decrease locally. In living organisms, this local decrease in entropy is typically achieved by the consumption of energy (more precisely, free-energy 1 ) and the increase of entropy in the surrounding environment, mostly in the form of heat and mass ows from the organism towards the environment. In this process, the main quantity of interest is the freeenergy, provided it can be meaningfully de ned, because it expresses the maximum amount of energy available to do work [Kondepundi and Prigogine, 2015]. Hence, in an alternative perspective, this free-energy expenditure is used to produce work that can counteract locally the entropic forces that tend towards disorder. Combining theses two ideas, we are led to another seemingly contradictory thought: Free-energy consumption contributes to the increase of entropy but it can also be used to produce work that acts against the entropic forces. Therefore, the crucial concept is that living beings funnel their entropy increase towards their environment, while keeping an ordered internal state 2 .

In living beings, the energy expenditure required to maintain homeostasis (the state of steady internal conditions in living beings), can be related to metabolism. However, an energy-consuming metabolism is not the only characteristic of the living [Nelson, 2004]. The ability to reproduce is also widely considered as a hallmark of life. But energy consumption is also central to reproduction [START_REF] England | Statistical physics of self-replication[END_REF] and, more generally, to the key idea of evolution: natural selection (though the explicit connection, in this case, remains unclear).

In the context of a cell (which is the setting that will occupy us for the rest of this work), free energy is most commonly stored in the form of a chemical potential di erence between Adenosine Triphosphate (ATP) and Adenosine Diphosphate (ADP), often referred to as the energy currency of the cell [START_REF] Phillips | Physical biology of the cell[END_REF]. Indeed, ATP is a high-energy but stable molecule -it typically does not hydrolyse spontaneously in an aqueous solution [Westheimer, 1987] -that can be used to drive chemical reactions or as a source of work.

Therefore, it seems clear now that one of the characteristics of life is its ability to maintain itself away from thermodynamic equilibrium. However, our comprehension of non-equilibrium thermodynamics is far from complete. Mostly due to the pioneering work of L. E. Boltzmann and J. W. Gibbs, for over a century we have had a good understanding of the thermodynamics of systems at equilibrium [Gibbs, 1879]. Slowly, progress was made, by L. Onsager and I. Prigogine among others, regarding the behaviour of systems near equilibrium, within the framework of Linear Irreversible Thermodynamics (LIT) [Onsager, 1931, Glansdor andPrigogine, 1954]. More recently, a new area of research has emerged that focuses on the thermodynamic properties of P 3 small (molecular) stochastic systems far away from equilibrium [Jarzynski, 1997, Crooks, 1999, Seifert, 2005, Sekimoto, 2010], which has also enabled some progress in the understanding of chemical reaction networks at the macroscopic scale [START_REF] Rao | Nonequilibrium thermodynamics of chemical reaction networks: Wisdom from stochastic thermodynamics[END_REF]. As a result, systems inherently out of equilibrium, such as living beings, are starting to be understood at the microscopic scale; but there are still plenty of open questions at the macroscopic scale, due to, among other reasons, the spectacular complexity of living beings and the emergent behaviour that arises at larger length-scales [Anderson, 1972].

Nevertheless, order -a central concept to this thesis -can also exist at thermodynamic equilibrium, as is the case, at the atomic scale, in certain metallic alloys [START_REF] Ashcroft | Solid state physics. Holt, Rinehart and Winston[END_REF]; or, at the mesoscopic scale, in colloidal suspensions with an interplay between short-and longrange interactions [Campbell et al., 2005, Tarzia andConiglio, 2006]. In these cases, the system is considered to be in contact with a thermal reservoir, which allows for the energy exchange between the system and the environment, and enables the appearance of order by minimising the free energy of the system. Therefore, a balance between the energy of the system and its entropy appears, and order is produced when the system is capable of minimising its energy and lower its entropy while increasing that of the bath, ensuring that the overall entropy of the universe is still increasing.

Thus, there is not a single universal source of order in nature. While living beings seem to prefer the dissipation of free energy to maintain their ordered structures, there are other mechanisms by which order emerges without the need of dissipative structures. Then, several questions relevant to the topic of this thesis emerge: To which extent does organisation and morphogenesis within living beings require the consumption of free energy? Which tasks could be accomplished without dissipation? Far from answering these questions, we will content ourselves by describing the interplay that can appear when both sources of order act simultaneously. On the one hand, we consider passive mechanisms (no free energy is consumed), which, for the remainder of this work, are the interactions between solutes in a solution. On the other hand, we take into account active processes that do require free-energy dissipation. The types of active processes that we consider in this thesis are chemical reactions driven by a free-energy di erence between fuel and waste components (e.g., in the context of a biological cell, the chemical potential di erence between ATP and ADP).

Much of this thesis is hence devoted to the study of the rich phenomenology that can arise from the crosstalk between these two mechanisms. The system that we bear in mind throughout this work is that of a solution, which is not necessarily dilute, where interactions between the solutes may arise; and fuelled chemical reactions take place, driving the solution out of thermodynamic equilibrium. In contrast, much of the theoretical work in chemical reaction networks has considered only ideal solutions, which assumes that solutes are dilute and interactions among them, negligible. However, order can also emerge when these interactions are taken into account, as explained above, and the interplay between passive and active ordering mechanisms can give rise to novel self-organised structures. We believe that this point of view provides a more realistic, yet still oversimpli ed, description of the cytoplasm of a cell, given that it is often crowded with macromolecules [Roberts et al., 2002, Hö ing andFranosch, 2013] and can even segregate some P of its components in droplet-like condensates as a consequence of the interactions among solutes [START_REF] Brangwynne | Germline P granules are liquid droplets that localize by controlled dissolution/condensation[END_REF]. This type of solutions, hereafter termed as non-ideal solutions, are the cornerstone of this work.

Finally and on a more personal note, despite this introduction, which may give the impression of a well-planned research agenda undertaken during these 3 years and 3 months, most of this thesis is the result of ruthless opportunism. Rather than following a plan, the work behind the majority of the di erent chapters in this thesis was undertaken in a fairly anarchic way, with chance and curiosity playing a major role. While the topics treated here revolve around the same concepts and I have tried to construct a coherent line of thought throughout the di erent parts of this thesis, traces of this unplanned opportunism can be observed in the text that follows. P I:

I -

Thermodynamic equilibrium may be characterized by the minimum of the Helmholtz free energy de ned usually by: F = E -T S.

Are most types of "organisations" around us of this nature? It is enough to ask such a question to see that the answer is negative. Obviously in a town, in a living system, we have a quite di erent type of functional order.

I. R. Prigogine, Nobel Lecture, 1977 5

The advances in microscopy and quantitative measurements in cell biology has allowed to obtain huge amounts of quantitative data from experiments. This data, while intriguing and often spectacular, is not very informative by itself and it only attains its entire potential when placed within a theoretical framework, when confronted with hypotheses and ideas. Traditionally these frameworks where developed by biologists and were fundamentally qualitative or descriptive. However, as more quantitative data is produced, there seems to be a need for a robust mathematical framework that can make the most of such experimental output. In developing an appropriate quantitative theory, the ideas of theoretical physics, while not always directly applicable, can be of great help.

The rst part of this thesis is devoted to the mathematical modelling of two processes that control the spatio-temporal organisation in biological cells, which is established at the cellular scale by the collective behaviour of particles at the molecular scale. First, we describe mathematically a case of protein aggregation process whose dynamics are slowed down for large agglomerates due to the presence of intracellular obstacles. Our theoretical predictions are compared with experiments performed with an optogenetic protein. In the second half of this part, we model the bacterial nucleoid (that is, the part of a bacterial cell in which the chromosome is located) and explain its size and positioning. As mentioned in the Preface, our approach has two main ingredients: we consider the e ect of crowders and steric interactions and the consequences of non-equilibrium processes in the cell, such as, transcription and translation. In line with the idea pioneered by I. Prigogine [Prigogine, 1977, Kondepundi andPrigogine, 2015], these non-equilibrium processes are the source of order in the system and result in what he called dissipative structures.

Ideally, these two examples could be modelled starting from a common theoretical framework that applies to all dissipative structures. However, in practice, the particular details of each system make it a very challenging task and often it is best to explore particular models for each of the processes separately.

Nevertheless, the modelling of these two systems has a number of common features. One of them is the fact that their time evolution is derived from physical principles, such as, Fick's law of di usion or excluded-volume e ects. Speci cally, in both cases, the e ect of intracellular obstacles modulates the dynamics and localisation of the di erent objects in the system. In the cytoplasm of an eukaryotic cell, there seems to be a well de ned mesh size due to the cytoskeleton and other intracellular obstacles [START_REF] Luby-Phelps | Hindered di usion of inert tracer particles in the cytoplasm of mouse 3T3 cells[END_REF], Etoc et al., 2018], which in turn, as we will prove, regulates the size and growth of protein aggregates. In the bacterial nucleoid there is also a characteristic mesh size [START_REF] Xiang | Interconnecting solvent quality, transcription, and chromosome folding in Escherichia coli[END_REF], but in this case formed by the folded chromosome. This mesh size controls, among other things, the spatial location of translation within the cell, as messenger RNAs (mRNAs) with one or several translating ribosomes attached are excluded from the nucleoid to a great extent [START_REF] Sanamrad | Single-particle tracking reveals that free ribosomal subunits are not excluded from the Escherichia coli nucleoid[END_REF], probably due to their size being larger than the mesh size of the nucleoid.

As mentioned above, both systems are out of thermodynamic equilibrium, as is usual in living beings, which are characterised by, among other features, their continuous free-energy consumption. However, these two systems are out of thermodynamic equilibrium in di erent ways. While the protein aggregation system is out of equilibrium because their initial condition is not an equilibrium state, the bacterial nucleoid is out of equilbrium since there is a continuous expenditure of free energy to maintain a non-equilibrium steady-state. Therefore, in our rst example we study the relaxation dynamics towards equilibrium of protein aggregates and, in the second one, the characteristics of the non-equilibrium steady state. Finally, both processes control the spatial organisation on the cell but they do so in di erent ways: The dynamics of protein aggregation determine, as a function of time, the number and size of the aggregates, while in the bacterial cell the rate of transcription can control the non-equilibrium steady-state and, thus, the positioning of the nucleoid and the localisation of the di erent components of the transcriptional-translational machinery.

Given the di erences between the processes studied in this part of the thesis, the mathematical formalism used in each chapter is markedly di erent from the ones used in the other chapters. While in the second chapter (regarding the bacterial nucleoid) the analysis takes explicitly into account the types of interactions present in the system and it was mostly limited to its steady state, in the rst chapter (about protein aggregation) we have resorted to a simpli ed kinetic description of the system to study its time evolution in an analytically tractable way. These two approaches can be seen as two faces of the processes that living beings use to control their internal self-organisation: Time dependent processes, like growth or coalescence (e.g. [START_REF] Wang | Self-organized periodicity of protein clusters in growing bacteria[END_REF], Al Jord et al., 2021]), and patterned non-equilibrium steady-states ( [START_REF] Murray | Self-organization and positioning of bacterial protein clusters[END_REF] or [START_REF] Halatek | Rethinking pattern formation in reactiondi usion systems[END_REF]). Taken together, they enable the spatial organisation of biological cells with an impressive degree of complexity.

Chapter 1

Protein aggregation in cells

Protein aggregation is a process which spans multiple orders of magnitude both in time and space1 : From nucleation, when a couple of monomers of a given chemical species of interest bind together to initiate the process, to the formation of large clusters containing up to millions of monomers each [START_REF] Zidar | Characterisation of protein aggregation with the Smoluchowski coagulation approach for use in biopharmaceuticals[END_REF]. Protein aggregation is a widely studied topic in molecular and cell biology due to its connections with numerous diseases and disorders, such as, Alzheimer's and Huntington's disease [Selkoe, 2004]. In these diseases, large micron-sized aggregates appear in cells which are often associated with the loss of function of the proteins involved [START_REF] Ross | Protein aggregation and neurodegenerative disease[END_REF]. However, most importantly, protein aggregation has also been linked with the creation of toxic intermediate-size aggregates that interact inappropriately with functional cellular components [START_REF] Knowles | The amyloid state and its association with protein misfolding diseases[END_REF]. Furthermore, protein aggregates have been shown to exclude chromatin and disrupt the expression of certain genes, and the aggregate size appears to play an important role in this process [START_REF] Li | Real-time imaging of Huntingtin aggregates diverting target search and gene transcription[END_REF].

Thus, protein aggregation can be seen as a self-assembly process where monomers aggregate to form typically brillar structures, with important consequences for cell physiology. However, since toxic intermediates can be created in the process, the nal assembly of the system (a single large and inert aggregate, if we are considering irreversible aggregation) is not as important as its time evolution towards this nal state. Therefore, here we analyse a process where the kinetics are crucial for its understanding and the nal steady state is not as informative as the trajectory followed towards it. More precisely, we consider a process where aggregation is always thermodynamically favoured but can be kinetically arrested, due to the presence of intracellular obstacles.

In this chapter, we rst introduce the physico-chemical basis of the process of protein aggregation, describing the di erent kinetics that may arise in this type of systems and motivating our research. Then we proceed to brie y describe the experiments performed and the mathematical modelling of the process. Finally we compare the output of these two approaches and present the consequences of our research.

Physical origin of protein aggregation

Most proteins in the cell have a stable folded conformation which allows them to appear in a soluble form in the cytoplasm. However, certain proteins do not fold into a globular and soluble structure and, among them, a certain subset is prone to aggregation [START_REF] Knowles | The amyloid state and its association with protein misfolding diseases[END_REF].

Typically, a seed or nucleation event (which gathers a number of proteins to create an oligomer or nucleus of the aggregate) is needed to start the process, after which the protein oligomers grow and become large aggregates [START_REF] Narayanan | A rst order phase transition mechanism underlies protein aggregation in mammalian cells[END_REF]. This is due to the shape of the free energy of the oligomer: If the free energy of the oligomer as a function of its size in monomers (N ) has a single maximum at N = N * , then N * is the critical size of the nucleus and determines the smallest growth-competent unit [Šarić et al., 2016]. Therefore, there is a free-energy barrier to form a nucleus of size N * , after which the growth of the aggregate is thermodynamically favoured. While the nucleation step can be reversed (small oligomers may be created and afterwards dissolved, especially if the concentration of the protein of interest is not high enough), once large aggregates are formed we can normally consider the aggregation process to be irreversible to a great extent. It is in this later stage of the process that we will focus our e orts.

Given the importance of the kinetics in protein aggregation, a plethora of mathematical models have been developed, depending on the particular features of the aggregating system. For example, for brillar assemblies, good agreement was found between the experimental data and a mathematical model considering the elongation and fragmentation of brils in only one dimension [START_REF] Knowles | An analytical solution to the kinetics of breakable lament assembly[END_REF]. Here, however, we will not focus on the shape of the aggregate but on how its overall size and di usivity a ects its growth, as discussed below.

Di usion-limited irreversible aggregation in the cell

We have already established that, due to energetic considerations, once a large enough protein oligomer has been nucleated, it will continue to grow spontaneously and irreversibly, provided the right conditions are met, i.e. the solution is supersaturated [START_REF] Narayanan | A rst order phase transition mechanism underlies protein aggregation in mammalian cells[END_REF]. However, the kinetics of this growth are still very informative and important from the biological point of view as they determine the speed at which these aggregates will grow, their size and their number. Indeed, if the kinetics are slow enough, aggregation may never happen in a biologically relevant timescale, while if they are fast enough, few of the toxic intermediate-sized oligomers will exist (as they will soon become large and potentially inert [START_REF] Knowles | The amyloid state and its association with protein misfolding diseases[END_REF]). As stressed above, it is not only the nal state that matters but also the time evolution of the process.

There are two main factors controlling the timescale over which protein aggregates grow. One of them is the binding a nity, that is, provided two particles meet, how energetically favoured is their binding and, thus, how likely they are to bind. The other factor is the di usion coe cient, which determines the frequency with which particles bump into each other. If the limiting factor for the aggregation process is the binding a nity of the particles, we talk about reaction-limited aggregation, while if the limiting factor is the di usivity of the particles then the process is termed di usion-limited aggregation.

In this chapter, we consider, from the theoretical viewpoint, a di usion-limited irreversible 1.1. P aggregation process. As a result, the speed at which particles di use and meet entirely regulates the speed of the aggregation process. This theoretical approach is complemented by an experimental study (carried out by Aléria Duperray-Susini and Mathieu Coppey) which uses an optogenetic protein that upon blue light exposure oligomerises (CRY2olig). Given that, under frequent blue light exposure, the fragmentation of a cluster is strongly suppressed, one may view this process as an irreversible, out-of-equilibrium aggregation, which leads to the formation of increasingly large clusters. In other words, irreversibility stems from a negligible fragmentation or dissolution rate of the aggregates, at least in the timescale of the experiments [START_REF] Taslimi | An optimized optogenetic clustering tool for probing protein interaction and function[END_REF]. The role of di usion is also crucial, because it sets the speed at which the aggregation process unfolds: aggregation processes in the cellular cytoplasm are in uenced by the presence of physical obstacles which alter the di usion dynamics of growing aggregates within the cell [Wojcieszyn et al., 1981, Hö ing andFranosch, 2013]. In particular, in Ref. [START_REF] Etoc | Non-speci c interactions govern cytosolic di usion of nanosized objects in mammalian cells[END_REF], it was shown that, for quasi-spherical nanoparticles in HeLa cells, the di usivity drops by two or three orders of magnitude as the diameter of the nanoparticle is increased from 50 to 75 nm, due to steric interactions between the particle and the cytosolic meshwork of the cell. Particles above this thresholdwhich is sometimes referred to as the pore size of the cytoplasm -experience almost no di usion [START_REF] Luby-Phelps | Hindered di usion of inert tracer particles in the cytoplasm of mouse 3T3 cells[END_REF].

Thus, in view of these facts, we developed a theoretical and analytically tractable model with as few as two free parameters to mathematically describe an aggregation process where agglomerates of di erent sizes have very di erent di usivities. In collaboration with experiments, the model allows to validate the picture of two di erent timescales controlling large aggregate number and size. Overall, the successful integration and feedback between theory and experiments led to quantitative predictions and their veri cation, like the critical size above which an aggregate can be thought of as "trapped" in the cytoplasm.

For irreversible aggregation processes, the only possible steady-state is the one where all proteins form a single cluster. However, this is rarely the case in biological cells, as the cytoplasm typically exhibits multiple protein or enzymes clusters scattered all over its volume [START_REF] An | Reversible compartmentalization of de novo purine biosynthetic complexes in living cells[END_REF], Narayanan et al., 2019]. Therefore, the physical mechanisms which set the cluster number and size still remain a subject of investigation [START_REF] Buchner | Clustering and optimal arrangement of enzymes in reaction-di usion systems[END_REF], Castellana et al., 2014], and could have important biological consequences.

In the rest of the chapter, we address the problem of di usion-limited irreversible aggregation processes in cells, both theoretically and experimentally. We found that two timescales control aggregation processes: One related to the fast di usion of small clusters, and the other one, slower, potentially related to the hindered di usion due to the presence of intracellular obstacles. By matching theoretical predictions with experimental results, we estimated the threshold between these two timescales. However, on long timescales a di usive-like movement can be observed even for large aggregates, which was related to uctuations that stem from the incoherent e ect of a network of active forces in the cell [START_REF] Guo | Probing the stochastic, motor-driven properties of the cytoplasm using force spectrum microscopy[END_REF], such as the rearrangement of the cytoskeleton and endomembranes. Therefore, large aggregates, with a radius comparable or larger than this threshold, can be thought of as being strongly con ned, and subject to a dramatic hindrance in the di usivity which may have an important e ect on the aggregation dynamics. Consequently, the slow timescale above is likely to be determined by active uctuations that a ect the dynamics of intracellular objects larger than the typical pore size of the cytoplasm. Overall, our results shed light into the interplay between aggregation processes, and the dynamics of the crowded environment in the cell cytoplasm.

Experimental setup

The experimental system under study is an optogenetic protein CRY2olig, which oligomerises upon blue light [START_REF] Park | Optogenetic protein clustering through uorescent protein tagging and extension of CRY2[END_REF], fused to the uorescent tag mCherry, that is transfected into RPE1 cells (retina pigmented epithelium 1, mammalian cells). An important feature of this optogenetic protein is the persistence of its oligomerised state even in the dark, with a half life of around 23 mins [START_REF] Taslimi | An optimized optogenetic clustering tool for probing protein interaction and function[END_REF].

Twenty-four hours after transfection, cells are exposed to blue light: this blue-light exposure can be regarded as an out-of-equilibrium process, which triggers protein oligomerisation. The period at which cells are exposed to blue light is 120s, which is signi cantly lower than the half life of the oligomerised state in the absence of light stimulation (23 mins), allowing us to consider the aggregation process as irreversible. The dynamics of these protein clusters are then followed for one hour as shown in Figure 1.1 A: From the images we extracted initial protein concentrations, nal concentration of clusters and its size. For details regarding the experimental setup and image analysis, please see Appendix A.

The number of monomers cannot be determined directly from these images as there is a constant relating arbitrary intensity units and the actual concentration in each pixel. An estimate of this constant was obtained by imaging droplets with known concentration of mCherry (the uorescent tag used in the experiments) and comparing both images.

Qualitatively, we were able to distinguish two di erent dynamical regimes. The rst one is a regime characterised by rapid di usion and aggregation of small oligomers, see Fig. 1.1 A, which takes place right after the blue light is switched on, and lasts for a time lapse of the order 1.3. M 13 of minutes which is short compared to the imaging time of 1 hr. The second regime is characterised by larger clusters which exhibit slower di usion or almost no di usion, resulting in a slower aggregation process. Presumably, the lack of di usion is produced by the obstacles in the cytoplasm, which trap large aggregates. Hence, any di usion or movement of these large clusters could be widely attributed to movements and remodelling of the cytoplasmic structure [START_REF] Guo | Probing the stochastic, motor-driven properties of the cytoplasm using force spectrum microscopy[END_REF], Etoc et al., 2018]. These features are summarised in Figure 1.1 B. A similar drop in the di usivity between large protein aggregates and monomers was also observed in Ref. [START_REF] Li | Real-time imaging of Huntingtin aggregates diverting target search and gene transcription[END_REF].

The analysis of these images allowed us to obtain the cluster-size distribution, the cluster concentration and the mean cluster size as functions of the initial protein concentration, which are shown in Fig. 1.2.

Modelling the dynamical aggregation process

The theoretical basis of irreversible aggregation processes was introduced by von Smoluchowski over a century ago, and it can be summarised into his well-known equation for the coagulation kinetics [START_REF] Krapivsky | A kinetic view of statistical physics[END_REF]:

dc i (t) dt = 1 2 j+k=i k j,k c j (t)c k (t) -c i (t) k k i,k c k (t), (1.1)
where c i (t) refers to the concentration of clusters with i monomers, and k i,j is the aggregation rate constant (also known as aggregation kernel) between two clusters of mass i and j. The rate of aggregation k j,k c j (t)c k (t) is given by the law of mass action. Upon an appropriate choice of the aggregation kernel k i,j , Eq. (1.1) can adequately describe di usion-limited aggregation processes.

However, the kernel typically does not take account of the e ect of obstacles or pores, such as the ones found in the cytoplasm of a cell [START_REF] Luby-Phelps | Hindered di usion of inert tracer particles in the cytoplasm of mouse 3T3 cells[END_REF], Etoc et al., 2018]. To take account of this e ect, we leverage the insights from the experiments (previously described in Section 1.2) to build a kernel k i,j based on the separation of the two timescales involved in the aggregation process: On the one hand, there is a fast aggregation timescale (characterised by a rate α), involving monomers and small clusters that di use rapidly, and which ultimately leads to the formation of larger agglomerates. On the other hand, there is a slow aggregation timescale (with characteristic rate β) that comprises aggregates larger than the pore size of the cytoplasm. The threshold between these two timescales is the time, t 0 , beyond which all clusters are larger than the pore size of the cytoplasm. We denote by m * the cluster mass at which the agglomerate attains the size of the pore of the cytoplasm and barely di uses, see Fig. 1.1.

The objective of our model is to examine the e ect of a sharp drop in di usivity with particle size and, therefore, we neglect other hydrodynamic e ects, such as size-dependent di usivity. These assumptions are supported by the conclusions of Ref. [START_REF] Etoc | Non-speci c interactions govern cytosolic di usion of nanosized objects in mammalian cells[END_REF] where it was found that most of the drop in di usivity with particle size takes place in a narrow window of size and other variations in di usivity are small in comparison to this drop. In addition, this simple timescale-separation assumption allows to keep the complexity of the model low, while capturing the essence of the dynamics.

C

1. P

Fast-aggregation timescale

For the fast-aggregation timescale, we choose the following kernel for the Smoluchowski coagulation equation:

k j,k = α [θ(m * -j) + θ(m * -k)] , (1.2)
where α is a fast aggregation rate constant, and is assumed to be much larger than the rate constant associated with the slow aggregation timescale β, and θ(x) is the Heaviside step function. Since, within this timescale, clusters with mass larger than m * do not di use, they do not contribute to the aggregation kernel (1.2). In addition, we will assume that, by the time all clusters are of mass m * or larger, the contribution of the slow process to the clustering dynamics is negligible, therefore e ectively decoupling the timescales involved in the problem.

For the kernel of the type of Eq. (1.2), the following change of variables is known to simplify the Smoluchowski equation [Leyvraz, 2003]:

ϕ i (t) = c i (t)/N (t), dτ = N (t)dt (1.3)
where N (t) = i c i (t) and ϕ i (t) represents for the fraction of clusters that is of size i, which yields the following form for the Smoluchowski equations:

dϕ i (τ ) dτ = α j+k=i ϕ j (τ )ϕ k (τ )θ(m * -j) -ϕ i (τ )θ(m * -i) .
(1.4)

An important feature of Eq. (1.4) is its recursive structure, i.e., the equation for ϕ i only depends on ϕ j for j < i. One can prove inductively that the solution for the mobile clusters (ϕ i for i ≤ m * ) is given by:

ϕ i (τ ) = i k=1 (-1) k-1 e -kατ i -1 k -1 , (1.5)
where we have assumed that only monomers are present at t = 0, which mathematically it means that ϕ i (0) = δ i,1 (also known as monodisperse initial condition).

Proof of the solution (1.5)

Given that Eq. (1.4) is a recursive equation for ϕ i , in what follows we will attempt an inductive proof of the solution for any ϕ i with i ≤ m * , for which the Heaviside step function is equal to one. We will show that, for i ≤ m * , if the ansatz (1.5) holds for ϕ 1 , • • • , ϕ i-1 , then it holds for ϕ i as well. To achieve this, we insert the ansatz (1.5) in Eq. (1.4), where we evaluate

i-1 j=1 ϕ j (τ )ϕ i-j (τ ) = i-1 j=1 j n 1 =1 i-j n 2 =1 (-1) n 1 +n 2 -2 e -(n 1 +n 2 )τ α j -1 n 1 -1 i -j -1 n 2 -1 = i-1 j=1 i s=2 s-1 n 2 =1 (-1) s-2 e -sτ α j -1 s -n 2 -1 i -j -1 n 2 -1 , (1.6) 1.3. M
where, in the last equality, we have made the change of variable s = n 1 + n 2 . Now we can apply Vandermonde's identity:

m + n r = r k=0 m k n r -k , (1.7) which yields i-1 j=1 ϕ j (τ )ϕ i-j (τ ) = i-1 j=1 i s=2 (-1) s-2 e -sτ α i -2 s -2 =(i -1) i s=2 (-1) s-2 e -sτ α i -2 s -2 . (1.8)
Equation (1.4) for ϕ i , i ≤ m * now reads, assuming the ansatz (1.5) for ϕ j , j < i,

dϕ i (τ ) dτ = α(i -1) i s=2 (-1) s-2 e -αsτ i -2 s -2 -αϕ i (τ ) (1.9)
which can be rewritten as follows

d(ϕ i (τ )e ατ ) dτ = α(i -1) i s=2 (-1) s-2 e -α(s-1)τ i -2 s -2 (1.10)
and solved by direct integration along with the monodisperse initial conditions (which make the constant from the integration vanish), yielding .11) thus proving that if the ansatz (1.5) holds for ϕ 1 , • • • , ϕ i-1 , then it holds for ϕ i as well. It can easily be checked that the ansatz (1.5) is a solution of Eq. (1.4) for i = 1, 2 or 3 and, thus, it will hold for all integer i such that 0 < i ≤ m * .

ϕ i (τ ) = i s=1 (-1) s-1 e -sατ i -1 s -1 , ( 1
Alternatively, Eq. (1.11) can be recast into the form:

ϕ i (τ ) = (1 -e -ατ ) i e ατ -1 (1.12)
by the binomial theorem.

Evolution of the number of clusters

The number of clusters as a function of our rescaled time τ , during the fast aggregation timescale, can be obtained rst by summing Eq. (1.1) over all i > 0 with the kernel (1.2):

i dc i (t) dt = α 1 2 i,j+k=i c j (t)c k (t) [θ(m * -j) + θ(m * -k)] -α i,k c i (t)c k (t) [θ(m * -i) + θ(m * -k)] = α i,j+k=i c j (t)c k (t)θ(m * -j) -2α i c i (t) m * k c k (t) = -αN (t) m * j=1 c j (t).
(1.13) C 1. P and then by using the change of variables (1.3), which yields the equation for the cluster concentration in the rescaled time τ dN (τ )

dτ = -αN (τ ) m * i=1 ϕ i (τ ). (1.14)
Using the so-called hockey-stick identity:

m * n = m * -n k=0 k + n -1 n -1 , (1.15)
after inserting the solution (1.5) into Eq. (1.14) we obtain

dN (τ ) dτ = -αN (τ ) m * i=1 i j=1 (-1) j-1 e -jατ i -1 j -1 = -αN (τ ) m * j=1 m * -j î=0 (-1) j-1 e -jατ î + j -1 j -1 = -αN (τ ) m * j=1 (-1) j-1 e -jατ m * j , (1.16)
where, in the second equality, the following change of variable has been performed: i = î + j to alter the order of summation and, in the third equality, the hockey-stick identity has been used. The solution to Eq. (1.16) is

N (τ ) = N (0) exp    m * k=1 (-1) k-1 k e -kυα m * k τ 0    , (1.17)
where the brackets denote the di erence between their argument evaluated at υ = τ , and at υ = 0.

In the limit τ → ∞, which is equivalent to the concentration of clusters after the fast aggregation timescale has nished, the result is:

N (τ → ∞) = N (0) exp - m * k=1 (-1) k-1 k m * k , (1.18)
which can be rewritten in a way that makes it more explicit the dependency on m * :

N (τ → ∞) = N (0) m * exp - m * k=1 (-1) k-1 k m * k + log m * . (1.19)
We can identify

m * k=1 (-1) k-1 k m * k (1.20)
as the m * -th harmonic number, H m * , which diverge logarithmically as m * → ∞. The Euler-Mascheroni constant, γ, is the di erence between the n-th harmonic number and the logarithm of n in the limit where n → ∞: Hence,Eq. (1.19), in the large m * limit, can be recast as

γ = lim n→∞ n k=1 (-1) k-1 k n k -log n . (1.21) 1.3. M
N (τ → ∞) m * →∞ = N (0) m * e -γ . (1.22)
It is worth noticing that in the large-time limit the results only depend on the critical cluster mass m * . Thus, Eq. (1.22) is independent of the fast aggregation constant α, which implies that, for the model to be consistent, we only need that the timescale separation hypothesis α β is satis ed. By using Eq. (1.22), we obtain the fraction, Ñ , of clusters left in the system after the fast aggregation process has nished

Ñ (t 0 ) = N (t 0 ) N (0) = e -γ m * .
(1.23)

In Eq. (1.23) we have used that m * is large, which is justi ed by the fact that we expect m * to be of the order of 10 2 or 10 3 , and even for m * = 10 2 the error due to this approximation is less than 1%. We will now make use of Eq. (1.23) as the initial condition of the slow-timescale dynamics.

Slow-aggregation timescale

Minutes after the start of the aggregation process, the fast aggregation process is nished (around time t ∼ t 0 ). Based on the experiments, we assumed t 0 to be small compared to the nal time of the experiment t f = 1 hr, which allows us to neglect it and assume that the slow-aggregation timescale lasts for 1 hr (and not t f -t 0 ).

To describe the slow-aggregation regime t > t 0 , we write a Smoluchowski coagulation equation with k i,j = β, where β is the slow-aggregation rate. The solution in this case was rst given by von Smoluchowski [START_REF] Krapivsky | A kinetic view of statistical physics[END_REF]. Starting from the Smoluchowski coagulation equation (1.1), with the kernel k j,k = β:

dc i (t) dt = β 2 j+k=i c j (t)c k (t) -βc i (t) k c k (t), (1.24)
we sum over all i and rewrite the constraint on the summation j + k = i with a Kronecker delta, yielding

i dc i (t) dt = β 2 i,j,k δ j+k,i c j (t)c k (t) - i βc i (t) k c k (t) = β 2 j,k c j (t)c k (t) - i βc i (t) k c k (t), (1.25)
and with the de nition N (t) = i c i (t) we obtain

dN (t) dt = - β 2 N (t) 2 .
(1.26)

The solution to this equation is (B) we plot the cluster concentration and mean cluster size as functions of the protein density ρ measured in each cell, respectively. Grey dots correspond to results for individual cells, black squares to the average over 20 cells, and error bars to standard deviations. Blue lines correspond to least-square ts of the theoretical expressions for the cluster concentration and cluster size, that is, Eqs. (1.28) and (1.29). The orange dashed lines correspond to the predictions for a passive cytoskeleton, i.e., in the absence of an active dynamics (β = 0). (C) Cluster-size probability density function (PDF) for di erent protein concentrations: Low (below 300 monomers/µm 3 ), medium (between 300 and 600 monomers/µm 3 ) and high (above 600 monomers/µm 3 ). The black dashed line corresponds to the estimate of the cluster-mass threshold between timescales, m * , inferred from (A) and (B).

N (t) = 1 βt/2 + C , ( 1
where, imposing the initial condition N (0) = ρ Ñ (t 0 ) with ρ is the initial density of monomers, C takes the value (ρ Ñ (t 0 )) -1 . Integrating Eq. (1.26) from t = 0 to t f we obtain

N (t f ) = ρ Ñ (t 0 ) ρ Ñ (t 0 )t f β/2 + 1
.

(1.28)

Equation (1.28) has two unknown parameters: β and Ñ (t 0 ) which can be estimated from the experimental data. In particular, Ñ (t 0 ) can be obtained from relation (1.23), which implies that it is fully determined by the size at which the aggregate attains a size comparable to the pore of the cytoplasm (m * ). Furthermore, using the relationship ρ = N (t f ) m(t f ) (where m(t f ) stands for the mean cluster mass at time t f ), one can estimate the mean cluster mass at the end of the experiment as a function of the density of protein:

m(t f ) = ρ t f β/2 + Ñ (t 0 ) -1 .
(1.29)

In this chapter, equations (1.28) and (1.29), for the slow timescale, and Eq. (1.23) for the fast one, constitute our main theoretical results.

Comparison with experiments and conclusions

In order to test the model, we t Eqs. (1.28) and (1.29) to the experimental data for the mean cluster mass and cluster density at the end of the experiment (for details on the tting procedure, see Appendix A.3). Results are shown in Figure 1.2 A and B. The t yields β = 9.6 hr -1 µm 3

1.4. C and Ñ (t 0 ) = 1.4 × 10 -3 . Using Eq. (1.23), one obtains a value of m * = 390 monomers for the mass threshold above which clusters are expected to be trapped in the cytoplasm.

The experimental data also allowed us to quantify the cluster-size distribution, shown in Fig. 1.2 C. Our estimate of m * is close to the peak of the cluster-size distribution: This result is consistent with the assumptions in the model, where the aggregation process is slowed down for clusters of mass above m * .

We can also assess the consistency of our result with other experimental data by estimating the pore size of the cytoplasm from the prediction for m * . To this end, we can utilise the framework of di usion-limited cluster aggregation (DLCA) [START_REF] Kolb | Scaling of kinetically growing clusters[END_REF]. This theory assumes that clusters di use freely and bind to each other as soon as they come to contact. If the bonds created by each binding event are rigid and maintain their shape the resulting structure would be a very sparse fractal aggregate. However, in many cases this may not be true and the bonds may rearrange to make a more compact structure.

Taking into account this rearrangement, it has been suggested from numerical simulations that the fractal dimension of the DLCA clusters is approximately d f = 2.18 [START_REF] Meakin | The e ects of restructuring on the geometry of clusters formed by di usion-limited, ballistic, and reaction-limited cluster-cluster aggregation[END_REF]

, i.e., R m * r 0 d f ∼ m * (1.30)
where R m * is the radius of an aggregate of mass m * and r 0 = 2.5 nm is the radius of an individual CRY2olig monomer, that is, the average size of a protein containing ∼ 500 residues, which allows us to obtain a characteristic radius of the aggregate R m * . This radius stands for the typical size of the aggregates, and it does not strictly represent the radius of an aggregate, nor implies that the aggregate has spherical shape. Unlike other aggregating proteins that tend to form brils [START_REF] Qiang | Structural variation in amyloid-β brils from Alzheimer's disease clinical subtypes[END_REF], it is unclear the geometry of the structures formed by CRY2olig, which is why our estimate does not include any structural information about the aggregate other than the assumption that monomers bind to each other as soon as they come into contact.

Therefore, we obtain for the radius of an aggregate of mass m * : R m * ∼ 39 nm. It should be noted that these calculations are correct up to a constant that we cannot determine: Nevertheless, our estimates are coherent with the threshold found in Ref. [START_REF] Etoc | Non-speci c interactions govern cytosolic di usion of nanosized objects in mammalian cells[END_REF], where they set the threshold between di using and non-di using particles to be between 25 and 37.5 nm.

Finally, this estimate of the size, together with the good agreement between functional form predicted by the theory and the experimental data, allows us to validate the picture of two di erent timescales setting the number and size of large protein aggregates, to a great extent. It is possible to re ne this model by choosing more complex kernels k i,j for the Smoluchowski coagulation equation in order to explain the cluster-size distribution, but this simple model with only two free parameters (m * and β) already captures much of the dynamics and the empirical results.

Conclusion and outlook

In this chapter, we have studied di usion-limited aggregation of CRY2olig protein in mammalian cells, combining an experimental and a theoretical approach.

Our main result is the identi cation of two di erent timescales in the aggregation process: C 1. P On the one hand, there exists a short timescale where small clusters can freely di use and aggregate, leading to the formation of larger agglomerates. On the other hand, later on, large aggregates barely di use or do so very slowly. Based on previous work [START_REF] Luby-Phelps | Hindered di usion of inert tracer particles in the cytoplasm of mouse 3T3 cells[END_REF], Guo et al., 2014, Etoc et al., 2018], this e ect could be largely due to con nement of the aggregates within the cytoskeleton and other cytosolic obstacles. As a result, large clusters cannot di use nor aggregate, unless the con ning obstacles move or rearrange on a longer timescale. The predicted threshold between the two timescales corresponds to cluster sizes of ∼ 400 monomers, or ∼ 39 nm of radius, which roughly corresponds to the cytosolic pore size [START_REF] Etoc | Non-speci c interactions govern cytosolic di usion of nanosized objects in mammalian cells[END_REF].

Our model yields a quantitative estimate of the aggregation rate, β, relative to the long time scale: This rate would characterise the incoherent dynamics of an intracellular network of active forces, such as molecular motors [START_REF] Guo | Probing the stochastic, motor-driven properties of the cytoplasm using force spectrum microscopy[END_REF], which could thus be regarded as an active stirring of the aggregates.

In addition, our analysis demonstrates that clustering of CRY2olig in mammalian cells is markedly di erent from aggregation in a passive material with a xed pore size, where the dynamics of the aggregation process would halt as soon as the aggregates' size reaches the pore size. This comparison was made in Fig. 1.2 A and B, where the orange dashed lines represent the predictions for a passive material with the same pore size as that of the cells in our experiment (β = 0), while solid blue lines represent our model prediction, which includes the active stirring of clusters.

The ideas developed in this study can be generalised to a variety of biological systems that reach a steady state driven by non-equilibrium processes, such as synthesis, degradation, tra c or recycling of proteins [START_REF] Turner | Nonequilibrium raftlike membrane domains under continuous recycling[END_REF]. In addition, the mechanisms identi ed here could be extended to the kinetics of other intracellular phenomena, such as liquid-liquid phase separation [Garcia-Jove [START_REF] Navarro | RNA is a critical element for the sizing and the composition of phase-separated RNA-protein condensates[END_REF]. Indeed, systems under binodal phase separation might exhibit as well two timescales in their coarsening dynamics. The fast timescale rate, α, would represent the di usion-limited coalescence of droplets in the early kinetics. On the other hand, as droplets grow and di usion slows down, the main driving force of coarsening would presumably be Ostwald ripening, whose details could be taken into account by a parameter, or function, equivalent to the slow aggregation rate, β. Given that there is a free-energetic cost for a droplet to deform around a network of obstacles [START_REF] Shin | Liquid nuclear condensates mechanically sense and restructure the genome[END_REF], the e ect of obstacles in di usion would become important only for droplets with a characteristic radius R m * or larger. Thus, we expect the value of R m * for this case to be similar to the one predicted by our analysis.

Here, we have focused our analysis on the kinetics of protein aggregation and found that these kinetics too can in uence patterns of spatial localisation in cells, such as the size of protein aggregates. Therefore, this mechanism, by which the dynamics are arrested to a large extent, also plays an important role in setting the spatial organisation at the scale of the cell, at least in cases where the processes are slow compared to the timescales of interest in cells and, thus, a steadystate is never reached. In the following chapter we will study the converse case: that of a system that does reach a non-equilibrium steady state which organises the intracellular medium.

Chapter 2

Intracellular bacterial organisation: the nucleoid

Living systems show a high degree of organisation at multiple scales, from the molecular one to the macroscopic scales of organisms and ecosystems1 . The study of self-organisation and patterning in living beings has attracted much attention from a variety of scienti c communities, from molecular and developmental biology to physics, chemistry and applied mathematics. A seminal contribution from the latter eld is the morphogenetic mechanism based on reactiondi usion systems proposed by A. M. Turing in 1952 [START_REF] Turing | The chemical basis of morphogenesis[END_REF]. Nevertheless, in recent years, many other self-organisation mechanisms have been proposed such as mechano-chemical feedbacks [START_REF] Boocock | Theory of mechanochemical patterning and optimal migration in cell monolayers[END_REF], growth-related patterning [START_REF] Wang | Self-organized periodicity of protein clusters in growing bacteria[END_REF] and liquid-liquid phase separation (LLPS) [START_REF] Brangwynne | Germline P granules are liquid droplets that localize by controlled dissolution/condensation[END_REF]. In this chapter we will focus on the spatial organisation of bacteria, which usually lack membrane-enclosed organelles yet achieve a striking degree of organisation.

A notable example of spatial organisation in bacterial cells is that of their DNA: Despite the absence of a nuclear membrane, in many bacteria such as Escherichia coli (E. coli), the chromosome is not randomly spread throughout the intracellular space, but is markedly localised [Lewis, 2004, Bakshi et al., 2012], and forms a compact structure.

The E. coli chromosome is a circular double stranded DNA molecule of around 4.6 × 10 6 base pairs, ∼1.5 mm of contour length, with a persistence length of ∼ 50nm. Therefore, if the chromosome behaved as a random coil it would occupy a volume of ∼ 500 µm 3 [START_REF] Verma | Architecture of the Escherichia coli nucleoid[END_REF]. However, the typical E. coli cell is a few micrometers in length (a volume of ∼ 2µm 3 ), which means that the chromosome has to be highly compacted to t into the cell. In fact, not only does the chromosome t into an E. coli cell but it does not even ll it completely, as it localises to a particular region of the cell -the nucleoid [Lewis, 2004, Surovtsev and[START_REF] Surovtsev | Subcellular organization: A critical feature of bacterial cell replication[END_REF].

The localisation and degree of con nement of the nucleoid varies with growth rate and among bacterial species [START_REF] Gray | Nucleoid size scaling and intracellular organization of translation across bacteria[END_REF]. This organisation and localisation of the chromosome has C 2. O been shown to play an important role in many biological processes, including transcription via the distribution of RNA polymerases [START_REF] Weng | Spatial organization of RNA polymerase and its relationship with transcription in Escherichia coli[END_REF], translation via the localisation of mRNAs and ribosomes [START_REF] Sanamrad | Single-particle tracking reveals that free ribosomal subunits are not excluded from the Escherichia coli nucleoid[END_REF], and the positioning and di usion of protein aggregates [START_REF] Coquel | Localization of protein aggregation in Escherichia coli is governed by di usion and nucleoid macromolecular crowding e ect[END_REF].

In this chapter, we aim at theoretically describing the E. coli nucleoid, its formation, compaction and localisation at the various stages of the cell cycle. In order to do that, we will construct a mechanism reminiscent of LLPS, based on the steric interactions present in the crowded cytoplasm and including the non-equilibrium e ects of transcription and translation, which demonstrates that localisation patterns on the cellular scale emerge spontaneously from microscopic features on a molecular scale. In particular, we show that the segregation of the nucleoid from mR-NAs and ribosomes can be understood in terms of excluded-volume e ects only, which implies that it is a consequence of equilibrium statistical physics, as in classical phase separations. In addition, our analysis shows that other dynamical features, such as nucleoid positioning, can be driven by the synthesis and degradation of mRNAs, making it a purely non-equilibrium feature. Finally, we compare these results with experimental data obtained from cells growing lamentously [Wu et al., 2019], either with or without chromosome replication, providing important physical and mechanistic insights.

Characteristics of the nucleoid: Compaction and localisation

Despite the importance of chromosome localisation, the physical causes and regulatory mechanisms of its con nement are still largely unknown. One of the causes of the compaction of the nucleoid could be the fact that the cytoplasm acts as a poor solvent for the chromosome [START_REF] Xiang | Interconnecting solvent quality, transcription, and chromosome folding in Escherichia coli[END_REF], but many other factors could also a ect nucleoid compaction, like nucleoidassociated proteins that modify the folding conformation of the chromosome [Dorman, 2013].

The prokaryotic nucleoid is not as structured as its eukaryotic counterpart. Nevertheless, there are still a wide variety of proteins and enzymes that regulate certain aspects of it. One of them is the histone-like nucleoid-structuring protein H-NS that drives the formation of loops between the sites to which it binds [START_REF] Dame | Chromosome organization in bacteria: mechanistic insights into genome structure and function[END_REF], e ectively bridging between di erent locations in the chromosome. While H-NS primarily bridges between chromosome regions and controls loops, other nucleoid-associated protein bend the DNA, such as HU or Fis (factor for inversion stimulation) and have also been related to compaction of certain regions of the nucleoid [START_REF] Remesh | Nucleoid remodeling during environmental adaptation is regulated by HU-dependent DNA bundling[END_REF].

Furthermore, Dps (a nucleoid-associated protein whose expression increases drastically in stationary phase E. coli cells) is also known to cause global nucleoid compaction by unespecically binding to the DNA in a cooperative manner [START_REF] Almiron | A novel DNA-binding protein with regulatory and protective roles in starved Escherichia coli[END_REF], Janissen et al., 2018]. In contrast, the protein MatP binds speci cally to a motif called matS localised speci cally in the Ter region of the chromosome, compacting and organising this part of the genome only [START_REF] Mercier | The matp/matS site-speci c system organizes the terminus region of the E. coli chromosome into a macrodomain[END_REF].

Finally, the supercoiling of the chromosome can be altered by the action of topoisomerases (enzymes that break and rejoin DNA strands) [START_REF] Dame | Chromosome organization in bacteria: mechanistic insights into genome structure and function[END_REF]. For example, the enzyme gyrase has been shown to control supercoiling of the DNA plectoneme, which has been linked to 2.2. M nucleoid compaction [START_REF] Stuger | DNA supercoiling by gyrase is linked to nucleoid compaction[END_REF].

However, there are also physical forces that compact the nucleoid. It is often suggested that crowding agents in the cytoplasm exert a force against the chromosome which is a crucial factor for the compaction of the nucleoid [START_REF] Cabrera | Active transcription of rRNA operons condenses the nucleoid in Escherichia coli: Examining the e ect of transcription on nucleoid structure in the absence of transertion[END_REF], Zhang et al., 2009, Yang et al., 2020]. More generally, when an E. coli cell is lysed, their intracellular components spread across the solution and, in particular, its chromosome was observed to expand abruptly [START_REF] Pelletier | Physical manipulation of the Escherichia coli chromosome reveals its soft nature[END_REF], implying that it is mostly repulsive forces who hold the together nucleoid within the cell. Overall, this suggest that repulsive interactions like steric or excluded-volume e ects are crucial for nucleoid compaction -an insight that we will leverage in the following for the theoretical modelling of the nucleoid.

Regarding the localisation of the nucleoid within the cell, various studies have shown that, while the nucleoid is located at the center of the cell before chromosome replication, during and after replication the daughter chromosomes move out of the centre [START_REF] Joshi | Escherichia coli sister chromosome separation includes an abrupt global transition with concomitant release of late-splitting intersister snaps[END_REF], typically localizing at 1/4 and 3/4 positions on the long cell axis [Wu et al., 2019]. Indeed, this positioning is very robust and occurs for a very wide range of cell lengths [Wu et al., 2019], but it can be perturbed -or even destroyed -by treating the cells with several drugs, such as kasugamycin, due to yet unknown reasons [START_REF] Xiang | Interconnecting solvent quality, transcription, and chromosome folding in Escherichia coli[END_REF]. Moreover, the general mechanism underpinning this localisation pattern remains unclear and has been hypothesised to require an active process [Joyeux, 2019].

Modelling the nucleoid

Previous theoretical e orts to explain the compaction and localisation of the nucleoid have been mostly based on Monte Carlo simulations [START_REF] Mondal | Entropybased mechanism of ribosome-nucleoid segregation in E. coli cells[END_REF], Joyeux, 2019]. It was found that excluded-volume e ects between DNA and polysomes (mRNAs bound to multiple ribosomes) may account for segregation of the nucleoid from the rest of the cytoplasm [START_REF] Mondal | Entropybased mechanism of ribosome-nucleoid segregation in E. coli cells[END_REF].

By contrast, in this section, we develop a statistical-physics description of the spatial localisation of the molecular components of the E. coli transcriptional-translational machinery (TTM) -composed of DNA, mRNAs, and ribosomes -and identify a set of physical mechanisms underlying their localisation patterns. Unlike previous studies, we leverage semi-analytical methods, i.e., the virial expansion, which allow us to tackle the complexity of the system, and reduce it to a set of computationally tractable reaction-di usion equations.

We thus describe the dynamics of the E. coli TTM by means of a minimal non-equilibrium statistical physics model. The model includes steric interactions only -mathematically described by hard-core interacting potentials -because our aim is to provide a minimal framework that can capture the phenomenology of the bacterial nucleoid. Other interactions may be taken into account, such as weak van der Waals forces among macromolecules, the e ects of nucleoidassociated proteins or e ective attractive interactions between ribosomes and DNA, due to ribosomes which bind nascent mRNAs linked to DNA. However, there is evidence suggesting that excluded-volume e ects are the key players in the formation and morphology of the nucleoid [START_REF] Yang | The e ects of polydisperse crowders on the compaction of the Escherichia coli nucleoid[END_REF], Zhang et al., 2009], and that interactions between the chromosome and the cytoplasmic components are essentially repulsive [START_REF] Pelletier | Physical manipulation of the Escherichia coli chromosome reveals its soft nature[END_REF]. Given that in-C 2. O teractions are repulsive, to a rst approximation, they can be modelled by hard-core potentials [START_REF] Andersen | Relationship between the hard-sphere uid and uids with realistic repulsive forces[END_REF] and if there were any attractive interactions, provided these are not dominant, they could e ectively be taken into account by modifying the steric-interaction parameters (radius or length of particles). Indeed, as shown in the rest of this chapter, steric interactions alone can account for many distinct aspects of nucleoid behavior.

Apart from the inter-particle potentials, we also consider inherently non-equilibrium e ects, namely, transcription and translation. In a simpli ed approach, transcription is modelled by a rate of synthesis of mRNAs linear with the concentration of DNA. Transcription is described by the binding and unbinding of ribosomes to bare mRNAs or polysomes. Finally, we also include the degradation of mRNAs in the model as a linear term in the dynamical equations, see Eqs.

(2.1-2.3).

By observing that E. coli cells have an approximately cylindrical shape and symmetry, we reduce the three-dimensional cytoplasm to a single dimension along the long cell axis (see Fig. 2.1 A and B) and describe the TTM in terms of the one-dimensional concentrations of DNA segments, mRNAs, and ribosomes. Implicitly, this assumes that the cell is homogeneous in the radial direction, a hypothesis whose consequences will be assessed in the discussion section (Section 2.5). Within this homogeneity approximation, we denote by c DNA (x, t) the concentration of DNA plectoneme segments at position x along the long cell axis and time t, by ρ n (x, t) that of polysomes composed of an mRNA and n ribosomes, and by c F (x, t) that of freely di using ribosomes, see Fig. 2.1 A and B. We then consider the reaction-di usion equations for these concentrations, where we incorporate the di usion currents and the chemical reactions, i.e., ribosome-mRNA binding and unbinding, mRNA synthesis and degradation:

∂ t c DNA (x, t) = -∂ x J DNA (x, t),
(2.1)

∂ t ρ n (x, t) = -∂ x J n (x, t) -k on c F (x, t)ρ n (x, t) -k o n ρ n (x, t) + k on c F (x, t)ρ n-1 (x, t) + k o (n + 1)ρ n+1 (x, t) + α c DNA (x, t)δ n,0 -βρ n (x, t), (2.2) ∂ t c F (x, t) = -∂ x J F (x, t) -k on c F (x, t) n ρ n (x, t) + k o n nρ n (x, t) + β n nρ n (x, t).
(2.3)

In Eqs. (2.1-2.3), J DNA , J n , and J F denote the particle currents (derived in Subsections 2.2.1 and 2.2.3), k on and k o the rate constants for ribosome binding and unbinding due to completion of translation, respectively, α the rate at which mRNAs are created locally by transcription, and β the mRNA degradation rate. We assume no-ux (Neumann) boundary conditions for all the components of the TTM at x = ± where is the cell half-length.

Regarding the steric interactions, as shown in Fig. 2.1 A and B, we consider ribosomes as spheres of radius R and, because mRNAs and polysomes with n ribosomes are roughly globular polymer coils, we also approximate them as spheres of radius R 0 and R n , respectively. Because the E. coli DNA has a branched, plectonemic structure with a well-de ned persistence length and transverse radius [Odijk, 2000], we consider the chromosome as a set of cylindrical segments, where the length of each segment (L, the Kuhn length) corresponds to twice the persistence length [START_REF] Rubinstein | Polymer physics[END_REF]. For the sake of computational tractability, we treat the DNA segments as disconnected, as shown in Fig. 2.1 B.

Equilibrium free energy of an interacting inhomogeneous gas

Our model for the nucleoid assumes that all particles interact through hard-core potentials (also known as excluded-volume interactions) and the expression for the particle currents J must be obtained for such system. In order to do that we rst need to compute the free energy of the system at thermodynamic equilibrium.

We use the virial expansion, rst developed by Onnes [Onnes, 1902] over a century ago, to compute the free energy of a gas of interacting particles. Here we limit ourselves to state the result for hard-sphere potentials to second and third order in the expansion, and refer the interested reader to classical statistical-physics textbooks, e.g., [Huang, 1987, Pathria, 1996], for a complete explanation of the procedure, or to Appendix B.1 for a brief summary.

Free energy of a binary mixture of hard spheres

We rst limit ourselves to a binary mixture of hard spheres and later we will include the e ect of DNA cylinders. We consider N A and N B particles of species A and B, which are hard spheres with radii R A , R B and di usion coe cients D A , D B , respectively, con ned in a volume V. The hard-sphere potential V ij between the ith and jth particle is zero if the distance between particles is larger than the sum of their radii, and in nity otherwise. Then, the free energy of the system, up to second order in the virial expansion (third viral coe cients) is

F = -k B T log Z = k B T a=A,B N a log N a Λ 3 a V -k B T log 1 - a=A,B N a (N a -1) 2 V B (2) aa - N A N B V B (2) AB - a=A,B N a (N a -1)(N a -2) 6 V 2 B (3) aaa - a =b N a (N a -1)N b 2 V 2 B (3) aab , (2.4) C 2. O
where Z is the partition function as de ned in equilibrium statistical physics, k B is the Boltzmann constant, Λ a is the thermal de Broglie wavelength of species a and the ith virial coe cients B (i) ab are given by:

B (2) ab = 4π 3 (R a + R b ) 3 , (2.5) B (3) abc = 16π 2 9 R 3 b R 3 c + 3R a R 2 b R 2 c (R b + R c ) + R 3 a (R b + R c ) 3 + 3R 2 a R b R c (R 2 b + 3R b R c + R 2 c ) . (2.6)
By assuming that the volume is large, we can approximate the free energy as

F k B T a N a log N a Λ 3 a V + a=A,B N a (N a -1) 2 V B (2) aa + N A N B V B (2) AB + a=A,B N a (N a -1)(N a -2) 6 V 2 B (3) aaa + a =b N a (N a -1)N b 2 V 2 B (3) aab , (2.7)
where we have used the Taylor expansion of the logarithm.

Now we consider an in nitesimal distance dx, in which there are dN a molecules of each species, with one-dimensional concentration c a (x) = dN a /dx. The volume of each of these in nitesimal slices of the system is dV = σ dx, where σ is the cross section of the system. Then,

F 0 = - dF 0 - k B T dx a=A,B c a (x) log c a (x) + log Λ 3 a σ + a,b=A,B B (2) ab [c a (x)] 2 2σ + a,b,c=A,B B (3) abc c a (x)c b (x)c c (x) 6σ 2 , (2.8)
where we have approximated N a -1 with N a , which yields substantial simpli cation of the functional form of the free energy. In what follows we drop the term log Λ 3 a /σ because it only changes the chemical potential by a constant and will not alter the results. As a result, the logarithms with dimensional arguments remain.

The quantity dF 0 /dx in Eq. (2.8) is the free-energy density of the system assuming that the concentrations are uniform along the x axis, where this condition is denoted by the subscript '0'. If we assume that the local free energy density of an inhomogeneous system f = dF/dx, is a function of the uniform free-energy density and of the derivatives of the concentration, i.e., f = f (f 0 , ∇c a , ∇ 2 c a ...), then we can expand it around f 0 , considering the concentration and its derivatives as independent variables, as follows [START_REF] Cahn | Free energy of a nonuniform system. I. Interfacial free energy[END_REF]:

f = f 0 + a γ a d 2 c a dx 2 + k B T 2σ a,b κ ab dc a dx dc b dx , (2.9)
where κ ab are the Cahn-Hilliard coe cients, which account for spatial inhomogeneities in the concentrations and set the free energy cost of concentration gradients in the system. The second 2.2. M term in the right-hand side of Eq. (2.9) does not contribute to the total free energy of the system: when spatially integrated, this term vanishes because of the Neumann (no ux) boundary conditions. Then, the total free energy is

F = - dF = - dF 0 + k B T 2 -a,b κ ab dc a dx dc b dx dx, (2.10)
and the chemical potential for, e.g., species A reads

µ A (x) = δF δc A (x) = k B T 1 + log c A (x) + a=A,B c a (x)B (2) a,A + 1 2 a,b=A,B B (3) a,b,A c a (x)c b (x) - a=A,B κ a,A d 2 c a dx 2 .
(2.11)

For future convenience, we de ne ν a as the non-ideal contribution (which stems from steric interactions) to the chemical potential of each species, e.g., species A:

ν A (x) = a=A,B c a (x)B
(2)

a,A + 1 2 a,b=A,B B (3) a,b,A c a (x)c b (x) - a=A,B κ a,A d 2 c a dx 2 (2.12)

Free energy of the full TTM

We will now apply the ideas discussed above to the model for the TTM. In particular, we will use the virial expansion discussed for the binary mixture of hard spheres to obtain an expression for the free energy of the TTM. However, there is an important di erence that needs to be taken into account: not all the particles are spheres, because the DNA segments are considered to be cylinders of length L and radius ρ.

In addition, a DNA segment does not interact with another DNA segment in the same way in which it would interact with a polysome or ribosome. In fact, two overlapping DNA plectonemes may be nested into each other, as discussed in Ref. [START_REF] Mondal | Entropybased mechanism of ribosome-nucleoid segregation in E. coli cells[END_REF]. To model this nesting, while we use the radius ρ to describe overlaps between a DNA cylinder and ribosomes or mR-NAs in the virial expansion, we use a smaller, e ective radius ρ < ρ for overlaps between two DNA cylinders. We base the value of ρ on the hard-sphere model for DNA used for numerical simulations in Ref. [START_REF] Mondal | Entropybased mechanism of ribosome-nucleoid segregation in E. coli cells[END_REF], where each plectoneme segment is represented as a sequence of four bond beads and two node beads, and all beads have radius ρ. In the simulations, whenever a DNA segment collides with a particle which is not a DNA segment, none of the beads are allowed to overlap with the particle. On the other hand, whenever two DNA segments collide, node beads cannot overlap with each other, but bond beads can, according to the picture above. Given that the two node beads are located at the vertices which connect segments, each node bead contributes half of its volume to each plectoneme segment. The volume that a DNA segment excludes to other DNA segments, which we denote by πρ 2 L, is thus the volume of one node bead, Hard-core interactions between DNA segments (red cylinders) and ribosomes or polysomes (black spheres). DNA segments interact with each other through a cylinder of radius ρ , while they interact with polysomes and ribosomes with a radius ρ > ρ , satisfying Eq. (2.13).

i.e., one fth of the volume πρ 2 L that the segment excludes to particles other than DNA. As a result, we obtain the relation

ρ = ρ/ √ 5 (2.13)
between ρ and ρ . See Fig. 2.2 for a sketch of this interactions.

Moreover, ρ being an e ective parameter, it also includes e ects from the nucleoid-associated proteins that control chromosome folding at the molecular scale, as discussed in Section 2.1. Because nucleoid size can be a ected by varying the concentration and functionality of these proteins and enzymes, in theory it should be possible to modulate the value of ρ to e ectively account for these biochemical perturbations, at least in a phenomenological manner. Provided repulsive interactions are still dominant, even if there are molecular or enzymatic e ects in the compaction of the nucleoid, viewing the interactions within the cell as purely steric may help to obtain a more holistic view of nucleoid organisation, since we can integrate many di erent phenomena into a single parameter. Moreover, as long as the interactions are mostly repulsive, approximating these by a hard-core potential is theoretically well-grounded [START_REF] Andersen | Relationship between the hard-sphere uid and uids with realistic repulsive forces[END_REF].

Given the shapes of the particles discussed above, the functional form of the contributions of steric e ects to the virial expansion B (i) ab for ribosomes and mRNAs remains the same as for hard spheres, but the interaction of other species with DNA, and of DNA with itself, is di erent as we are considering DNA to be a set of disjoint cylinders. The second virial coe cient for two cylinders is [START_REF] Herold | Virial coe cients of anisotropic hard solids of revolution: The detailed in uence of the particle geometry[END_REF]

B

(2)

DNA DNA = 2πρ Lρ + 1 2 (L + ρ )(L + πρ ) , (2.14)
and the virial coe cient between one cylinder and a sphere with radius R n can be computed by performing integrals of the form (B.4) and yields

B

(2)

DNA n = Lπ(R n + ρ) 2 + 2πR n ρ 2 + πρR n 2 + 2 3 R 2 n , (2.15)
where the subscript 'DNA' stands for a DNA cylinder, and n for a sphere of radius R n .

Before we present the expressions for the third virial coe cients for cylinders, let us de ne

u ij = I(q i ∩ q j ), (2.16) 2.2. M 29
where q i denotes the degrees of freedom which specify the position and orientation of a particle, i.e., the position of its center of mass and the Euler angles de ning the particle orientation. In addition, q i ∩ q j in Eq. (2.16) stands for the condition that particles i and j overlap, i.e., their hard-core potential is nonzero, and the indicator function I is equal to one if the condition in its argument is satis ed, and zero otherwise.

The third virial coe cients for interactions that involve cylinders are given by the following integral expressions:

B (3) D n m = dr n D dr m D u z D m u z D n u nm ,
(2.17)

B ( 3
) D D n = 1 8π 2 dr D D dr D n dr D n d(cos θ D )d(cos θ D )d(cos θ n )dφ D dφ D dφ n u z D D u z D n u D n , (2.18) B ( 3
) DD D = 1 8π 2 dr DD dr D D dr DD d(cos θ D )d(cos θ D )d(cos θ D )dφ D dφ D dφ D u z DD u z DD u D D , (2.19)
where the subscripts label di erent cylinders, 'D' is a shorthand for DNA, n and m label the spheres, and vectors r ij denote the relative position between the centers of mass of particles i and j. In addition, the superscript z means that the axis of cylinder D is parallel to the z axis, so as to leverage spherical symmetry, and θ, φ are the polar and azimuthal angles, respectively. While some simpli cations of those integrals are possible, there is no known analytical form for these virial coe cients [Straley, 1973], and we obtained them by numerical integration. Nevertheless, this numerical integration is much less computationally costly than the simulations performed in Refs. [START_REF] Mondal | Entropybased mechanism of ribosome-nucleoid segregation in E. coli cells[END_REF], Joyeux, 2019], as it only needs to be performed once for xed parameters.

Then, the total free energy of the system is:

F k B T = - dx a=F,n,DNA B (2) DNA,a σ c DNA (x)c a (x) + a,b=DNA,F,n B (3) DNA,a,b 2σ 2 c a (x)c b (x)c DNA (x) + - dx dF 0 k B T + a,b=DNA,F,n κ ab 2 dc a dx dc b dx , (2.20)
where the sums run over the chemical species denoted by F, DNA and all polysome species, which we denote by 'n' in the sums. Note that F 0 has the same structure as Eq. (2.8): in fact, its form does not change because it involves only species of spherical particles. The di erence between F 0 in Eq. (2.20) and F 0 in Eq. (2.8) is in the summation indices, which now span over all polysome species and free ribosomes, which are all modelled as spheres.

The free energy of the system should include an entropic term related to the DNA conformation, but, since this term is small compared with the steric interactions, we neglect it for simplicity (a more detailed explanation is given in Section 2.2.2). The steric interactions have been written explicitly for the DNA segments, and the Cahn-Hilliard terms are analogous to those in Eq. (2.10), except for the numerical values of the coe cients κ ab , which depend on the particle geometry.

In principle, the Cahn-Hilliard coe cients κ ab in Eq. (2.20) can be computed by leveraging hard-core interactions as shown in Ref. [START_REF] Ilker | Phase separation and nucleation in mixtures of particles with di erent temperatures[END_REF]. However, for our purposes it is enough to observe that the Cahn-Hilliard terms re ect the cost of concentration gradients related Because κ ab has the dimension of the cube of a length, it must be proportional to a product of the linear sizes of the particles of species a and b. In addition, the relation of the Cahn-Hilliard coe cients to di erentials of concentrations over in nitesimal length scales indicates that κ ab is physically related to short rather than long length scales. We thus assume that κ ab equals the minimum between the volume of species a and that of species b.

DNA free energy

When deriving the particle currents, the quantity of interest is the free energy of the particles. For independent spherical particles (ribosomes and polysomes), the entropic term in the free energy is included in the virial expansion. However, for the DNA plectoneme the situation is di erent: Due to the lack of connectivity between the DNA "cylinders" in the model, we do not obtain the correct entropic term from the virial expansion.

In our model of DNA composed by disconnected cylindrical segments, the energetic part of the free energy is given by the interactions between DNA segments, encoded by the virialexpansion terms, see Subsection 2.2.1. Here, we argue that the entropic part of the real free energy can be neglected, because it is much smaller than the virial terms.

The main contribution is related to the free-energetic cost of con ning an ideal polymer. This entropic cost is given by [Edwards andFreed, 1969, De Gennes, 1979]:

S -k B N DNA L 2 L 2 n , (2.21)
where L n is the typical lengthscale on which the polymer is con ned and N DNA is the number of DNA segments. In the case of the nucleoid, L n ∼ 1 µm and N DNA = 6 × 10 3 .

By inserting the numerical values of the parameters we obtain an entropic contribution to the free energy of the order of 10 3 k B T and a contribution from the virial terms of order

k B T N 2 DNA B (2) DNA,DNA 2V n ∼ 10 4 k B T, (2.22)
where V n ∼ L 3 n is the typical volume of a nucleoid. Therefore, there is a di erence of one order of magnitude between the entropic and the virial term and, for the sake of simplicity, in this work we neglect the entropic contribution to the free energy, and thus to the current of DNA segments.

Auxiliary entropy

In our analysis of Eqs. (2.1) to (2.3), we will rst determine the steady state of the system in the absence of reaction and non-equilibrium terms, by numerically minimising the total free energy (2.20). These pro les are then used as initial conditions to integrate forward in time the reaction-di usion Eqs. (2.1) to (2.3), which include both reaction and non-equilibrium terms, see Subsection 2.2.5 for details on the numerical methods to solve Eqs. (2.1) to (2.3). At the freeenergy minimum, the DNA concentration is nonzero in the nucleoid, while it vanishes outside the nucleoid. Given that these equilibrium pro les are entered as initial conditions in Eqs. (2.1) to (2.3), the vanishing concentration above causes numerical instabilities when these equations are numerically integrated forward in time, and can lead to negative concentrations in the nonequilibrium steady state [Shampine et al., 2005]. To overcome this issue, we included a small, additional entropic term in the free energy (2.20). Therefore, this additional terms is included solely for the stability of the numerical methods used. For details about this term see Appendix B.2, and for a comparison between the magnitude of the real free energy and the auxiliary one see

Particle currents

We can now work out the currents by considering Fick's law of di usion [Crank, 1979], namely,

J a (x) = - D a k B T c a (x)∂ x µ a , (2.23)
where D a is the di usion constant, which is given by Einstein-Smoluchowski-Sutherland relation [START_REF] Dill | Molecular driving forces: statistical thermodynamics in biology, chemistry, physics, and nanoscience[END_REF]]

D a = ζ a k B T, (2.24)
and ζ a is the mobility of species a.

In Eq. ( 2.23) the chemical potential µ a is obtained from the derivative of the equilibrium free energy, see Eq. (2.11). We are thus constructing a hydrodynamic theory where we only consider variations in the slow modes, e.g. concentration di erences along the cell, implicitly assuming that uctuations in the fast modes are quickly thermalised and equilibrated, thus locally obeying equilibrium statistical physics with a well de ned local temperature, pressure and chemical potential. This is known as the assumption of local equilibrium [Kondepundi and Prigogine, 2015] and is the starting point of the Linear Irreversible Thermodynamics framework used to derive these currents [START_REF] Groot | Non-equilibrium Thermodynamics[END_REF].

Particle current for hard spheres

By substituting the expression (2.11) for the chemical potential in Eq. (2.23), we obtain the current for hard sphere particles, e.g. free ribosomes:

J F (x) = -D F ∂ x c F (x) + c F (x) a=F,DNA,n B (2) a,F ∂ x c a (x) + c F (x) 2 a,b=F,DNA,n B (3) a,b,F [∂ x c a (x)c b (x) + c a (x)∂ x c b (x)] -c F (x) a=F,DNA,n κ a,F d 3 c a dx 3 (2.25)
Similarly, we obtain the current for the rest of the spherical species (mRNAs and polysomes).

Particle current for the full TTM

Proceeding along the lines of the case for hard spheres, we obtain the DNA current from Eq. (2.23), using the free energy (2.20) and the auxiliary free energy (B.15):

J DNA (x) = -c DNA (x)D DNA [∂ x ν DNA (x) + ∂ x µ aux (x)] , (2.26) C 2. O
where ν DNA is the excluded-volume term analogous to that in Eq. (2.12), which stems from the fact that the chemical potential of DNA segments does not have entropic contribution. The term µ aux is the derivative of the auxiliary free energy (B.13) with respect to the DNA concentration, c DNA (x), and its contribution to the current reads

∂ x µ aux (x) = K aux e -Dauxc DNA (x)/ c DNA (x) ∂ x c DNA (x) 1 c DNA (x) 1 -D aux c DNA (x) c DNA (x) + D aux log[σc DNA (x)] - D aux c DNA (x) 1 + log[2 σc DNA (x)] 1 -D aux c DNA (x) c DNA (x) .
(2.27)

Combining the results obtained for hard spheres, Eqs. (2.8) and (2.25), with those for the full TTM, Eqs. (2.20) and (2.26), we obtain the currents in Eqs. (2.1) to (2.3). Making use of Eqs. (2.25) and (2.26) and of the virial coe cients previously derived, the currents for DNA, ribosomes and polysomes are fully de ned.

For completeness, the currents for ribosomes and polysomes are

J F (x) = - D F k B T c F (x)∂ x µ F , (2.28) J n (x) = - D n k B T ρ n (x)∂ x µ n , (2.29)
respectively, where the chemical potentials are

µ F (x) = log c F (x) + a=F,DNA,n c a (x)B
(2)

a,F + 1 2 a,b=F,DNA,n B (3) a,b,F c a (x)c b (x) - a=F,DNA,n κ a,F d 2 c a dx 2 ,
(2.30) for free ribosomes, and

µ n (x) = log ρ n (x) + a=F,DNA,n c a (x)B (2) a,n + 1 2 a,b=F,DNA,n B (3) a,b,n c a (x)c b (x) - a=F,DNA,n κ a,n d 2 c a dx 2 ,
(2.31) for polysomes.

Model parameters

We x the model parameters from experiments as follows. First, we consider the parameters on a molecular scale: The radius and length of DNA cylinders are ρ = 10 nm and L = 200 nm [Odijk, 2000, Mondal et al., 2011], respectively, where L/2 is approximately the persistence length of a DNA plectoneme2 [START_REF] Cunha | Polymer-mediated compaction and internal dynamics of isolated Escherichia coli nucleoids[END_REF], Odijk, 2000].

We take the ribosome radius to be R = 10 nm [START_REF] Mondal | Entropybased mechanism of ribosome-nucleoid segregation in E. coli cells[END_REF], and the radius of a ribosome-free mRNA to be R 0 = 20 nm [START_REF] Kaczanowska | Ribosome biogenesis and the translation process in Escherichia coli[END_REF]. The radius R n of an mRNA loaded with n ribosomes is estimated as the sum of the volume of a bare mRNA and n times the volume of a ribosome, i.e., 4/3π(R 3 0 + nR 3 ) yielding R n = (R 3 0 + nR 3 ) 1/3 . We estimated the di usion constant of the di erent species as follows: D F = 0.4 µm 2 /s for ribosomes, and D n = 5 × 10 -2 µm 2 /s for bare mRNAs and polysomes [START_REF] Bakshi | Superresolution imaging of ribosomes and RNA polymerase in live Escherichia coli cells[END_REF], Sanamrad et al., 2014]. For DNA segments it is harder to obtain an estimate based on experimental data. We therefore estimate their di usion coe cient as follows: We assume that the drag coe cient on the nucleoid as a whole is low enough for the nucleoid to relax relatively rapidly to mechanical equilibrium (for details, see Appendix B.4) and, given that DNA segments have a linear dimension similar to that of polysomes, we assume that their di usion coe cients will also be similar. Therefore, we take D DNA = 10 -2 µm 2 /s, a value for the di usion coe cient of DNA segments which allows the nucleoid to react rapidly to perturbations, but is still not larger than the di usion coe cient of polysomes.

The parameters relative to the cellular scale are the total number of ribosomes per cell N F , the cell half-length , both of which will be varied, and the radius of the cellular cross section, which is held constant. Because a central aim is to compare to the experiments in [Wu et al., 2019], we are interested in values for a doubling time of ∼ 2 hr as in that study. We thus interpolated experimental data points for di erent growth rates, to obtain the parameter values for the desired growth rate (see Appendix Sections B.5 and B.7) and obtained a total number N F ∼ 7300 ribosomes, a cross-sectional radius R cell ≈ 0.4 µm, and a cell half-length ∼ 0.9 µm for a reference cell. In addition, the total mRNA concentration for the reference cell was xed at ρ tot = n ρ n = 2400 µm -3 [START_REF] Bartholomäus | Bacteria di erently regulate mRNA abundance to speci cally respond to various stresses[END_REF]. The total number of DNA cylinders for the reference cell was taken to be N DNA ∼ 6700 segments [START_REF] Mondal | Entropybased mechanism of ribosome-nucleoid segregation in E. coli cells[END_REF]. When analysing di erent situations of biological interest, these parameters may be varied, particularly the overall concentrations of DNA segments, ribosomes, and mRNAs.

Finally, we set the reaction constants for the ribosome binding and unbinding to k on = 6 × 10 -4 µm/s, k off = 2.5 × 10 -2 /s [START_REF] Castellana | Spatial organization of bacterial transcription and translation[END_REF]. The mRNA degradation rate β = 3 × 10 -3 /s corresponds to an mRNA half life of ∼ 5 min [START_REF] Bernstein | Global analysis of Escherichia coli RNA degradosome function using DNA microarrays[END_REF], and the mRNA synthesis rate α is estimated from the global steady-state condition of Eq. (2.2), αN DNA = βN mRNA [START_REF] Castellana | Spatial organization of bacterial transcription and translation[END_REF], where N mRNA is the total number of mRNA molecules in the cell, i.e., ρ tot times the cell volume.

Numerically solving Eqs. (2.1) to (2.3)

Due to the complexity and nonlinearity of Eqs. (2.1) to (2.3), in order to solve them we need to resort to numerical methods. We will solve these equations in two stages: First we will minimise the free energy of the system to obtain an equilibrium steady state -which is not yet a solution to Eqs. (2.1) to (2.3) -and then we will integrate forward in time Eqs. (2.1) to (2.3) using the equilibrium steady state as initial condition, in order to obtain a non-equilibrium steady state. This technique for solving the equations re ects the two types of e ects we are considering. On the one hand, we have interactions that do not require any source of external work and, thus, in the absence of any other e ect, relax to thermodynamic equilibrium. This implies that their steady state can be found by minimising their free energy [Pippard, 1964]. On the other hand, C

2. O when integrating forward in time, we are adding intrinsically non-equilibrium e ects (since transcription and translation require the cell to consume free energy), which lead the system to a nonequilibrium steady state. Below we argue why a minimisation of the free energy yields a steady state. For further details on the numerical methods used and, in particular, for the integration in time of Eqs. (2.1) to (2.3) from an equilibrium steady state, see Appendix B.6.

Free-energy minimum and equilibrium steady state

The di usion equations for the mixture of hard spheres discussed in Section 2.2.3 (that is, in the absence of any non-equilibrium contribution) read

∂ t c a (x, t) = -∂ x J a (x, t), (2.32)
where a stands for each of the di usive species involved.

The steady state of Eq. (2.32) combined with no-ux boundary conditions, and with a constraint which xes the total number of particles to a given value (N a ) for each species yields the following set of relations: where F is now considered to be a functional of the concentration pro les c a (x), as discussed in Section 2.2.1, and is minimised with respect to the concentration pro les of all species a.

-∂ x J a (x, t) = 0, J a (± , t) = 0, (2.33) - dx c a (x, t) = N a . ( 2
The Lagrange function of the minimisation problem given by Eqs. (2.35) and (2.36) reads (2.37) where λ a are the Lagrange multipliers. First, the stationarity condition of L with respect to c a (x) is given by (2.38) where in the third line we used the de nition of the chemical potential, see the rst line of Eq. (2.11). By taking the derivative of Eq. (2.38) with respect to x and using Eq. (2.23), we obtain J a (x) = 0 for all x, which is equivalent to Eq. (2.33). Second, the stationarity condition of L with respect to λ a yields Eq. (2.34). As a result, the minimisation problem (2.35), (2.36) is equivalent to conditions (2.33) and (2.34).

L = F - a λ a - dx c a (x) -N a ,
0 = δL δc a (x) = δF δc a (x) -λ a = µ a (x) -λ a ,

N 35

Nucleoid formation and size

So far we have described the model we will be working with and we now turn to analyse the features of it. The rst characteristic of the model we study is whether it gives rise to a nucleoid segregated from the crowders (mRNA and ribosomes) and its size.

A system that relaxes to thermodynamic equilibrium at constant volume and in thermal contact with a constant temperature bath, minimises its free energy F = U -T S, where T is the temperature of the bath [Pippard, 1964]. In certain cases, the way a system minimises its freeenergy implies segregating some components of the mixture from others. This is the case of the mixture of DNA and crowders (mRNAs and ribosomes) which we are considering here, where, as it will be shown below, the DNA demixes from ribosomes and mRNAs. This phenomenon is known as phase separation and the region of space to which the DNA localises after the segregation can be identi ed with the bacterial nucleoid. Intuitively, in our model for the nucleoid, phase separation is a consequence of the tendency of the system to reduce its interaction energy (and, thus, its overall free energy), which is achieved by segregating plectonemic DNA from polysomes and ribosomes. However, it can also be viewed as steric interactions giving rise to depletion forces [START_REF] Asakura | On interaction between two bodies immersed in a solution of macromolecules[END_REF], which drive the phase separation of the nucleoid.

In this section we will give a brief introduction to phase-separation phenomena and analyse the equilibrium behaviour of the model. We will then numerically solve Eqs. (2.1) to (2.3), which include non-equilibrium transcription and translation, and compare these solutions with analytical estimates for the nucleoid size.

Segregation of DNA and crowders

Given a certain set of overall particle number constraints [e.g. the type of constraint introduced by Eq. (2.34)], the fact that the system relaxes to a free-energy minimum implies that, at that minimum, the Hessian matrix of the system

H i,j = ∂ 2 F ∂x i ∂x j , (2.39)
where x i are the degrees of freedom and F is the free energy of the system, has only positive eigenvalues. If we constrain the system to be homogeneous and nd that one or more eigenvalues of the Hessian matrix are negative, then the homogeneous con guration is unstable: focusing on the case of phase separation, this means that the system will demix into two or more distinct phases. This phenomenon is known as spinodal decomposition. In other words, whether spinodal decomposition occurs depends on the curvature of the free-energy function [De Gennes, 1979].

We now consider a simpli ed version of the free energy density (2.20) of the model:

f k B T =c F log c F + n ρ n log ρ n + B (2) F,F c 2 F 2 + n B (2) n,F c F ρ n + 1 2 n,n B (2) n,n ρ n ρ n + B (2) F,DNA c F c DNA + n B (2) n,DNA c DNA ρ n + 1 2 c 2 DNA B
(2) we can obtain the number of free ribosomes and polysomes in a homogeneous system from the kinetics of ribosome binding and unbinding, see Eqs. (2.1) to (2.3):

c F = c F,tot 1 + ρ tot kon k o , ρ 0 = ρ tot exp   - c F,tot k o kon + ρ tot   and ρ n = 1 n! ρ 0 k on c F k o n , (2.41)
where, in order to obtain these analytic results, we have assumed that there is no upper bound on the amount of ribosomes that can bind to an mRNA [START_REF] Castellana | Spatial organization of bacterial transcription and translation[END_REF], which is justi ed if the amount of ribosomes that can bind an mRNA is large enough. Substituting these relations in the free energy (2.40), we can obtain the spinodal region of the phase diagram by evaluating the eigenvalues of the Hessian matrix as functions of the concentrations and parameters of the system. If we choose the parameters as Section 2.2.4, but vary the total mRNA concentration ρ tot and the radius by which DNA segments interact with other DNA segments ρ (which are the two most relevant parameters for the phase behaviour of the system under consideration), we obtain the phase diagram shown in Fig. 2.3.

While the phase diagram Fig. 2.3 provides an estimate for the parameter region where phase separation takes place, the correspondence with the full model is not exact for a number of reasons. First, it is a simpli ed free energy where only second virial coe cients have been used for simplicity, as opposed to the full model, which includes also third-order virial coe cients. Second, the phase diagram is based on equilibrium thermodynamics, which need not apply to the nonequilibrium model de ned by Eqs. (2.1) to (2.3). In this sense, we have obtained the equilibrium phase diagram that is conceptually closest to the non-equilibrium system which we are modeling, because the phase diagram incorporates non-equilibrium e ects from ribosome binding and unbinding, as imposed by Eqs. (2.41). Nevertheless, our phase diagram cannot incorporate the e ects that may arise from mRNA synthesis and degradation, such as the centering and splitting of the nucleoid that, as it will be shown later, are consequence of transcription. Therefore, the phase diagram Fig. 2.3 is a qualitatively good approximation to the true non-equilibrium dynamics (e.g. it predicts phase separation for the parameter ranges explored numerically later in the text) but it is not exact.

N

Estimate for the spinodal decomposition within the nucleoid

Here, we will consider a further simpli ed version of the free energy of the system to estimate the values of the parameters needed to drive spinodal decomposition within the nucleoid, a quantity that will be useful later for the understanding of the splitting of the nucleoid.

We consider the simple case of only two species: DNA segments and bare mRNAs, with the following free energy

F = ρ 0 log ρ 0 + 1 2 B (2) 00 ρ 2 0 + B (2) 0DNA c DNA ρ 0 + 1 2 c 2 DNA B
(2)

DNADNA .

(2.42)

Setting all parameters as explained in Section 2.2.4, we can obtain the eigenvalues of the Hessian matrix with respect to these two variables (c DNA and ρ 0 ) and when any of the two eigenvalues becomes negative, spinodal decomposition occurs. Thus, we nd that the concentration of mRNA required at the center of the nucleoid to drive the splitting of the nucleoid is ρ * 0 2500 µm -1 (for the typical DNA-segment concentration within the nucleoid of c DNA 6000 µm -1 , see below and Fig. 2.4 A).

Again, we note that this is not an exact approach and is probably an overestimate of the required mRNA density at the center of the nucleoid. The fact that we only consider bare mRNA and not larger polysomes will cause some of this overestimate. Moreover, considering only second virial coe cients will also cause overestimation of the required mRNA concentration.

Nucleoid size in the presence of non-equilibrium processes

In order to better understand the compaction of the nucleoid, we solved the one-dimensional reaction-di usion Eqs. (2.1) to (2.3). To compare the predictions of our model with experimental data, in what follows we consider two scenarios for how the concentrations of the molecular species scale with cell length.

Filamentous growth

In lamentous growth, the total number of DNA segments, mRNAs, and ribosomes is proportional to the cell length. Given that one of the primary aims of our analysis is to compare results with the data of Ref. [Wu et al., 2019], which was obtained for slowly growing cells, the cell length is kept xed, which yields a substantial computational simpli cation. For each cell length, we rst determined the equilibrium steady state of the system by minimising the free energy (2.20), and then numerically integrated the reaction-di usion Eqs. (2.1) to (2.3) forward in time to reach a non-equilibrium steady state. The results are shown in Fig. 2.4 for di erent cell lengths, up to the cell length at which the nucleoid spontaneously splits into two lobes, and for di erent total mRNA densities. The cell length at which the nucleoid spontaneously splits into two lobes is easily identi able because, when the splitting takes place, the DNA concentration at midcell drops to a value close to zero, comparable to the DNA concentration at the cell poles.

While the non-equilibrium chemical reactions are responsible for the splitting of the nucleoid (as it will be argued in Section 2.4.1), the minimisation of the free energy can account for the existence of a phase-separated nucleoid in the cell due to the steric interactions that give rise to Cell length 2ℓ (µm) Nucleoid length and standard deviation as a function of cell length from Ref. [Wu et al., 2019] (red) and from the model (black). We found agreement up to cell lengths of around 9µm, after which the prediction of the model deviates from the experimental data.

L n (µm) αc DNA (µm s) -1 L n (µm) αc DNA (µm s) -1
depletion forces [START_REF] Asakura | On interaction between two bodies immersed in a solution of macromolecules[END_REF], see Fig B.4 for results obtained in the same way as in Fig. 2.4 but in the absence of non-equilibrium e ects. Importantly, neglecting the nonequilibrium e ects did not yield a substantial change in nucleoid size.

Furthermore, we observe that the con guration that minimises the free energy excludes ribosomes, mRNAs, and polysomes from the nucleoid to di erent degrees. For example, the free ribosome concentration is higher at the periphery than within the nucleoid by a factor of ∼ 3 while for a polysome with 3 bound ribosomes the ratio rises to ∼ 40, which is consistent with experimental observations that polysomes are much more excluded from the nucleoid than are ribosomes [START_REF] Sanamrad | Single-particle tracking reveals that free ribosomal subunits are not excluded from the Escherichia coli nucleoid[END_REF]. As a rule of thumb, we nd that the larger the macromolecular complex is, the more excluded it will be from the nucleoid.

The relation between nucleoid and cell length appears to be roughly linear up until the cell length at which the nucleoid begins to split in two, see Fig. 2.4 A and B. Moreover, as shown in Fig. 2.4 B, the higher the total mRNA density, the smaller the nucleoid, implying that a high mRNA density increases the osmotic pressure on the nucleoid, thus making it shrink. In Section 2.3.3 we will make a more formal statement of this observation and we will use it to estimate the nucleoid size analytically.

Single chromosome lamentous growth

So far we have analysed the scaling of nucleoid size with cell size by assuming that the number of DNA segments is proportional to cell length. We now study another case of biological interest, namely, that of a cell with a xed amount of DNA and varying cell size. This scenario was recently analysed in a dynamic imaging study of the E. coli chromosome [Wu et al., 2019], where the initiation of DNA replication and cell division were halted, yielding a single chromosome in 40 C

2. O a lamentously growing cell. We model this scenario by xing the number of DNA segments, but allowing the cell size to vary. In addition, the mRNA and ribosome number are no longer proportional to cell length: based on the data in [Kohram, 2021], we assume that the total concentrations of mRNAs and ribosomes decrease linearly with cell length (ρ tot ∝ a -b , where a and b are constants), approaching zero at 30 µm -see Appendix B.9 for details.

Results are shown in Fig. 2.5: the model again predicts a roughly linear scaling of the nucleoid size with respect to cell length, while the DNA-segment concentration decreases with cell size. This indicates that the decrease in DNA-segment concentration with cell size is balanced by the decrease of mRNA and ribosome concentrations, so as to keep nucleoid size a linear function of cell size. This can be seen clearly in Fig. 2.5 A where the concentrations of all components of the TTM decrease as the cell size increases. While the model prediction for nucleoid versus cell length agrees reasonably well with experiments [Wu et al., 2019] for cell lengths smaller than ∼ 10 µm (Fig. 2.5 B), there is a discrepancy for larger cells , see Section 2.5 for further discussion on these results.

Finally, in Fig. B.5, we show an alternative case in which the mRNA and ribosome concentrations are kept constant instead of decreasing with cell size. This yields a predicted nucleoid length very small compared to the experimental data, even for short cells, lending support to the results of Ref. [Kohram, 2021] and the scaling proposed above.

Analytical estimates of the nucleoid size

As seen above, the nucleoid segregates from the rest of the cytoplasm. The nucleoid size at steady state (provided that the nucleoid is single lobed) is set by the mechanical balance of osmotic pressures between the nucleoid and the peripheral cytoplasm. These pressures solely stem from the entropy and steric interactions of the components of the mixture, making the nucleoid size a consequence of equilibrium physics, which is in line with the conclusions of the previous section.

Motivated by the phase separation of the nucleoid and its exclusion of large components of the cytoplasm, in order to estimate the size of this phase-separated nucleoid, we will consider the cell to be divided by two movable walls into three compartments: A central one, the nucleoid, composed exclusively of DNA segments3 , and two lateral ones which include ribosomes and polysomes. To reach mechanical stability, the compartments may expand or contract, moving the walls to an equilibrium position where their osmotic pressures are balanced. The particles interact through steric interactions, described by the virial coe cients. Therefore, we consider the following free energy for the particles within the compartments:

F i = F ideal + k B T N 2 i B i 2V i , (2.43)
where i denotes the compartment and N i the particle number, B i the e ective virial coe cient that accounts for the steric interactions among the particles within the compartment, F ideal is the free energy of the ideal gas and V i the compartment volume. Essesntially, we are coarse graining all the second virial coe cients of the di erent particle types that appear in the full free energy (2.20) into a single one for the whole compartment B i . Then, the osmotic pressure exerted by the compartments is [De Gennes, 1979]

P i = - ∂F ∂V = k B T N i V i 1 + N i B i 2V i , (2.44)
where the rst term stems from the entropic pressure of an ideal gas while the second one,

N i B i /V i , comes from the steric interactions.
In what follows, we consider Eq. (2.44) in the nucleoid and in the pole compartments, and estimate the respective values of the virial coe cients in the nucleoid, B n , and at the cell poles, B p .

In the nucleoid, DNA-DNA interactions dominate, yielding a value of B n 6.4 × 10 -4 µm 3 . For the poles we provide an e ective value of the virial coe cient by assuming that all ribosomes are bound to mRNAs, and are equally distributed among them, that is, the poles are occupied by spheres all equal in size. Given that the ratio of ribosomes to mRNAs changes with the amount of mRNA in the cell, the virial coe cients depend on this last parameter. In the cases analysed in Fig. 2.4 B, we obtain the following values for the virial coe cients, using Eq. (2.5), for the corresponding mRNA concentrations: B p (1500 µm -1 ) = 4.3 × 10 -4 µm 3 , B p (2400 µm -1 ) = 3.7 × 10 -4 µm 3 , B p (3000 µm -1 ) = 3.5 × 10 -4 µm 3 .

(2.45)

By equating the pressures of the compartments, we obtain the equilibrium value for the volumes of each compartment. For lamentous growth, where the number of ribosomes, mRNAs, and DNA segments scales linearly with size (N i ∝ ), we obtain the solution for the nucleoid size V n = ϕV , where ϕ, the fraction of total volume occupied by the nucleoid, depends on the concentration of mRNAs in the cell: ϕ(1500 µm -1 ) = 0.7, ϕ(2400 µm -1 ) = 0.61, ϕ(3000 µm -1 ) = 0.56.

(2.46)

As shown in Fig. 2.4 B, these estimates (gray lines) are in good agreement with the numerical solution of the full reaction-di usion equations (red points).

Nevertheless, the cytoplasm of the cell may not be dilute enough to allow considering only two terms in the virial expansion. The third virial coe cient for the DNA-DNA interaction is B

(3) n 5 × 10 -8 µm 6 and for the polysomes: B p (1500 µm -1 ) = 9.4 × 10 -8 µm 6 , B p (2400 µm -1 ) = 6.8 × 10 -8 µm 6 , B p (3000 µm -1 ) = 6 × 10 -8 µm 6 .

(2.47) C

O

If we add the third virial coe cient, the fraction of total volume occupied by the nucleoid becomes:

ϕ(1500 µm -1 ) = 0.69, ϕ(2400 µm -1 ) = 0.60, ϕ(3000 µm -1 ) = 0.55, (2.48) values that are very close to those obtained with the second virial coe cient alone. Therefore, for simplicity, when comparing with the nucleoid size obtained by numerically solving Eqs. (2.1) to (2.3) we use the value obtained from second order virial coe cients only.

Overall, we nd that the inclusion of steric terms in both the nucleoid and mRNA/ribosome compartments makes the nucleoid swell compared to what its size would be with only entropic terms (ideal gas contribution). This is due to the nature of the nucleoid, a long relatively sti polymer with little entropy per segment compared to ribosomes and mRNAs. Moreover, in Fig. 2.4 B we quantitatively compared the estimates made in this section (grey lines) with the results of the numerical solution to the full model (red data points), which includes non-equilibrium e ects, and found good agreement between the two approaches; thus validating the idea of nucleoid size being a consequence of equilibrium physics alone and that the e ects of non-equilibrium transcription and translation in nucleoid size are negligible, provided the levels of mRNA and ribosomes are maintained. Therefore, nucleoid size is a consequence of steric interactions alone, which can be parametrised by an e ective osmotic pressure to which only large macromolecules contribute (mostly mRNAs and polysomes), since they are the only ones excluded from the nucleoid.

Consequences of non-equilibrium transcription and translation in intracellular spatial organisation

So far we have analysed the emergence, compaction and size of the nucleoid and concluded that they are a consequence of steric interactions, i.e., mostly driven by processes that do not require free-energy dissipation. In other words, they are passive mechanisms. In this section, we will explicitly address the consequences of non-equilibrium transcription and translation in nucleoid positioning, as a result of the particle currents that can be established at steady state.

Currents at steady state

When transcription and translation are taken into account new phenomena occur, particularly at steady state. As seen previously in this chapter, integrating numerically Eqs. (2.1) to (2.3) up to their non-equilibrium steady state yields a pro le like the one shown in Fig. 2.6. However, in Fig. 2.6, we stress the fact that there are di usion currents sustained at steady state due to the transcriptomic activity. Indeed, mRNAs are synthesised within the nucleoid with rate αc DNA (x), hence their spatial distribution has a local maximum at midcell. Then, mRNAs di use towards the poles (see ux arrows in Fig. 2.6), since it is free-energetically favoured. Finally, mRNAs reach the cell poles where, after a typical timescale ∼ β -1 = 5 min, they are degraded. Note that mRNAs within the nucleoid have no or few ribosomes attached while at the poles they are heavily loaded The plot depicts the concentration proles at the non-equilibrium steady state. The blue arrows in the bottom depict the steady-state ux of mRNAs, J ρtot (the arrow length is proportional to the mRNA current). At steady state, mRNAs are synthesised within the nucleoid, wherein they di use before escaping to the poles (and are ultimately degraded). with them, creating a circulating current for ribosomes too [START_REF] Castellana | Spatial organization of bacterial transcription and translation[END_REF]. Altogether, these processes lead to the steady-state distribution shown in Fig. 2.6. In Fig. 2.6 the arrows depict the ux of mRNAs, which is quanti ed by J ρtot (x) = n J n (x), with J n (x) given by Eq. (2.29).

These steady-state currents are at the heart of the spatial localisation of the bacterial nucleoid.

Nucleoid splitting

While the linear increase of nucleoid length with cell length is the result of equilibrium osmoticpressure balance, the splitting of the nucleoid is entirely due to non-equilibrium processes. In fact, for cells with ρ tot = 2400/µm and a half length of ∼ 4 µm or larger, the equilibrium steady state used as the initial condition for the reaction-di usion equations yields a nucleoid with a single lobe. By contrast, the nucleoid splits into two identical lobes positioned at 1/4 and 3/4 of the long cell axis when the reaction-di usion Eqs. (2.1) to (2.3) are integrated forward in time, see Fig. 2.7. Such 1/4 and 3/4 positioning of the daughter nucleoids has been ubiquitously observed in experiments [Wu et al., 2019] and is reproduced by the model with a high degree of accuracy, see Fig.

Scaling for the length at which the nucleoid splits

In what follows, we present a simple argument to explain the dependence of the length at which the nucleoid splits with respect to the underlying parameters, e.g., the mRNA synthesis rate. We take the nucleoid to be a region of length L nucl with homogeneous DNA-segment concentration which extends from x = -L nucl /2 to x = L nucl /2, with interfaces that are perfectly sharp. The mRNAs synthesised within the nucleoid di use until they reach the nucleoid boundaries and, because it is free-energetically favourable, they then escape the nucleoid and not return. As a result, the steadystate concentration of mRNAs within the nucleoid can be modeled by the following di usion equation with a uniform source term due to mRNA synthesis and absorbing boundary conditions, which represent mRNAs escaping from the nucleoid:

D n ∂ 2 ρ tot (x) ∂x 2 + αc DNA = 0, ρ tot ± L nucl 2 = 0, (2.49)
where ρ tot (x) is the total mRNA concentration at position x, D n the mRNA di usion constant (as de ned in Section 2.2.4), and αc DNA the rate of mRNA synthesis. The solution of the above equation is

ρ tot (x) = (L 2 nucl /4 -x 2 ) α c DNA /(2D n ),
whose local maximum within the nucleoid at x = 0 takes the value L 2 nucl α c DNA /(8D n ). We hypothesise that when the mRNA concentration at the center becomes larger than a given threshold, ρ * tot , spinodal decomposition takes place due to steric interactions between mRNAs and DNA, causing the nucleoid to split into two lobes. We thus expect ρ * tot to roughly correspond to the spinodal line of the phase diagram, but, given the non-equilibrium nature of the system mostly due to mRNA synthesis, it could di er from the equilibrium spinodal boundary. If we assume the threshold ρ * tot exactly corresponds to the spinodal threshold, then ρ * tot 2.5 × 10 3 µm -1 (estimated in Section 2.3.1). However, we note that an mRNA concentration of ∼ 10 3 µm -1 to drive the splitting of the nucleoid is in agreement with the observed behavior of the full model (see Fig. 2.7) and, therefore, we can justify a posteriori the validity of the approximations.

Whatever value ρ * tot takes (provided its dependency on α is negligible), this simple model predicts a scaling for the critical length L * nucl at which the nucleoid starts to divide of the form L * nucl ∝ (αc DNA ) -1/2 , obtained from equating the maximum of the mRNA concentration pro le to a xed value ρ * tot . To test the prediction of this simple model, we numerically obtained the length at which the nucleoid divides for di erent values of α, see the inset in Fig. 2.4 B, and found a good agreement with the proposed scaling.

Nevertheless, our model is too simple to give a mechanistic explanation for nucleoid splitting: no connectivity between DNA segments or modelling of DNA replication is included in the description. However, overall, the model does show that mRNA synthesis can control the localisation of daughter nucleoids, once the nucleoid has been split in two.

Nucleoid centring

As observed in Ref. [Wu et al., 2019], a single bacterial nucleoid has a strong tendency to localise at midcell for all cell sizes. Following the recent suggestion that the central positioning of the nucleoid is regulated by an active process [Joyeux, 2019], we investigated whether the non-2.4. C equilibrium process of mRNA production, di usion, ribosome binding, and mRNA degradation can account for nucleoid centring.

We consider the case of a nucleoid that, due to a uctuation, is not initially at the center of the cell, and test whether the non-equilibrium e ects in our model can push the nucleoid back to the cell center. To model this, we use the steady-state pro les obtained for lamentous growth, and shift the concentration pro les towards the right cell pole. The resulting con guration has a nucleoid displaced from the center, and equal mRNA and ribosome concentrations on both sides of the nucleoid. This concentration pro le is used as the initial condition for Eqs. (2.1) to (2.3), which we integrate forward in time in the presence of the non-equilibrium terms. As shown in Fig. 2.8, the nucleoid is centred at midcell after ∼ 30 min.

The physical origin of this centring is mRNA synthesis within the nucleoid: The nascent mRNAs di use in the nucleoid until they reach one of its boundaries and then escape, with an equal ux to the left and right of the nucleoid. If the nucleoid is not centred, the accumulating mRNAs occupy a greater fraction of the available volume on one side of the nucleoid and thus create a higher osmotic pressure on that side. The resulting pressure di erence ultimately drives the nucleoid back to the center of the cell.

The rate at which the nucleoid moves towards the cell center depends on both the pressure di erence due to mRNA accumulation, and on the e ective viscous drag experienced by the nucleoid as a whole. The e ective viscous drag, γ, is given by F nucl = γ v nucl , where F nucl is the force applied on the nucleoid (here, the osmotic-pressure di erence times the cross-section of the cell), and v nucl is the velocity of the nucleoid. If γ is low, the nucleoid responds fast to any force applied on it and thus quickly reaches a position of mechanical equilibrium, where the osmotic-pressure di erence vanishes. In this case, the centring process is only limited by the speed at which mR-NAs accumulate on either side of the nucleoid, which sets the pressure di erences. By assuming this limit, we can to establish a lower bound for the time needed by the nucleoid to center.

The kinetics obtained in this low-drag limit are shown in Fig. 2.8B (grey line), and they are given by an exponential relaxation with timescale β -1 , set by the rate of mRNA degradation (see below). As shown in the gure, the nucleoid centring obtained from the full model lags behind the lower bound, showing that there is a non-negligible contribution from drag on the nucleoid. As shown in the inset, both the lower bound and the result from the full model show an exponential relaxation of the nucleoid position for early times in the centring process.

Analytical estimates for the centring dynamics in the low-drag limit

The centring of the nucleoid can also be understood in terms of the simpli ed compartment model we used to estimate the nucleoid size. However, in order to explain the centring dynamics, we need to modify the model and assume that the mRNA synthesised in the nucleoid can di use out of the nucleoid to the lateral compartments. If the nucleoid is not centred in the cell and the synthesised mRNAs leave the nucleoid symmetrically to the left and right, then the mRNA density, and thus the osmotic pressure, will increase in the smaller polar compartment, thus pushing the nucleoid towards the center. As a result, the force that we need to consider is the di erence in pressure between the poles times the cross section σ of the cell F = σ(P L -P R ), where 'L' and 'R' denote the left and right pole, respectively. The dynamical equation for the position of the center of mass of the nucleoid, which we denote by x n , is:

dx n dt = D k B T σ(P L -P R ), (2.50)
where D is a di usion constant, not necessarily equal to the di usion constant of DNA segments.

In fact, D is an e ective di usion coe cient that includes collective e ects of DNA segments di using together and potentially other biological e ects.

We can provide a lower bound for the time it takes the nucleoid to center by assuming that the drag is small, i.e. the centring of the nucleoid due to a di erence in osmotic pressure is only limited by the synthesis of mRNA. In this limit, the nucleoid moves fast enough to prevent a pressure di erence between the poles, that is, a quasi-static approximation of Eq. (2.50): d t x n = 0, which implies P L = P R . Thus, the centring of the nucleoid is controlled by the rate at which the number of polysomes in the lateral compartments change. In both compartments, the pressure is set by the concentration of polysomes, whose number is set by the following di erential equation:

dρ tot i dt = αc DNA 2 -βρ tot i , (2.51)
where ρ tot i is the total mRNA concentration in compartment i and whose solution is

ρ tot i = α c DNA 2β + C i e -βt , (2.52)
where C i is a constant that is set by initial conditions. In the case of the initial condition of Fig. 2.8, C i takes the value C L +0.2 αc DNA /β and C R -0.2 α c DNA /β for the left and right compartment, respectively, as the initial position of the centre of mass of the nucleoid is located at +0.2/ and the amount of mRNA is directly proportional to the volume of each compartment. Since P L = P R we have ρ tot L /V L = ρ tot R /V R . Assuming that the nucleoid does not change size during this process, we obtain

V L = ρ tot L ρ tot L + ρ tot R (V -V n ), (2.53) 2.4. C
where V is the total volume of the cell, and V n the volume of the nucleoid. In the previous relation, the only term that is time-dependent is ρ tot L since ρ tot L + ρ tot R is constant in time.

Therefore, the position of the nucleoid is set by V L , which depends only on ρ tot L , which yields the exponential relaxation with timescale β -1 mentioned above. This lower bound on the time for centering is depicted in Fig. 2.8 B (grey line).

Validity of the one-dimensional model

In our analysis, we leveraged the cylindrical symmetry of the E. coli cell to reduce the number of dimensions of the model to one, which, in the description of the centring dynamics of the nucleoid, is crucial. The real system is three-dimensional and, as opposed to what has been assumed before, it need not be homogeneous in the radial direction. This becomes particularly important at midcell, where the nucleoid may not occupy the whole cellular cross section, potentially allowing the di usion of polysomes and ribosomes around the nucleoid: such di usion could have an important e ect on the centring dynamics, as it would result in a decrease of the osmotic-pressure di erence between the poles and, therefore, a decrease of the centring force. Actually, it would mean that the osmotic pressure di erence can vanish without requiring the nucleoid to move, thus disrupting the centring mechanism proposed above. However, the magnitude of this e ect depends on how fast polysomes can di use around the nucleoid. If such di usion around the nucleoid is slow (due to, e.g., obstacles in the narrow channel between the membrane and the nucleoid), the osmotic-pressure di erence, and thus the centring force, will be sustained.

We now proceed to estimate the di usive ux around the nucleoid in the presence of a concentration gradient and assess the e ect of this ux in the centring dynamics. We assume that the nucleoid is a cylinder of length L nucl and radius R -L/2, where R is the cell radius and L is the length of a DNA segment (we estimate the depletion zone between the membrane and the nucleoid to have a width equal to the persistence length of a DNA plectoneme as this is the lengthscale of the DNA-membrane interaction). The cross section available for di usion is

σ d = π[R 2 -(R -L/2) 2 ]
πRL (neglecting terms of order L 2 ). If we assume that the concentration of mRNAs at one pole is c and that it is zero at the opposite pole, then the concentration gradient is c/L nucl , yielding for the total ux across J the available cross section

J D n c L nucl πRL, (2.54)
where D n is the di usion coe cient of mRNAs and polysomes.

Then, the characteristic timescale τ over which these concentration gradients will disappear due to di usion is

τ ∼ cV p J = RL p L nucl D n L , (2.55)
where L p is the length of one of the polar regions (2 = 2L p + L nucl ) and V p = πR 2 L p is its volume. By substituting the values of these parameters given in Section 2.2.4, and assuming that for our reference cell L nucl 1 µm and L p 0.5 µm, we obtain τ ∼ 20 s.

The value obtained for the di usion timescale τ is lower than the typical timescale on which synthesis and degradation of mRNAs (β -1 ) would set the osmotic-pressure di erence at the C 2. O poles. Therefore, most of the osmotic-pressure di erence would be removed by di usion, and the resulting centering force exerted on the nucleoid would be small. However, when estimating τ , we assumed that polysomes di using through the narrow channel between the nucleoid and the membrane have same di usion coe cient as for polysomes di using at the poles. This assumption is likely to be inaccurate, because in the channel there are more crowders or obstacles than at the poles, including transertion e ects -that is, the simultaneous transcription, translation and insertion of proteins into the membrane that tethers the nucleoid to the membrane [START_REF] Gorle | DNA condensation in live E. coli provides evidence for transertion[END_REF]. Hence, the validity of the one-dimensional model depends on how fast diffusion of mRNAs and polysomes occurs through the narrow channel between the nucleoid and the membrane. For su ciently slow di usion in this region and a su ciently narrow channel, the corrections to the one-dimensional model (stemming from the nucleoid not occupying the whole cross section) are negligible.

Perturbing the nucleoid

Until now we have analysed how the model behaves naturally, in the absence of stresses or perturbations. However, we can also probe how the model would respond to certain types of perturbations and compare with experimental data. Usually, antibiotics are utilised to modify certain essential functions in E. coli and observe the changes in the size and morphology of the nucleoid.

Examples of these perturbations are a halt of transcription with rifampicin and an arrest in translation by kasugamycin. Both scenarios are analysed below.

Halt in transcription

It has been shown experimentally that when transcription in E. coli is halted, e.g. by treatment with rifampicin, the nucleoid expands [START_REF] Cabrera | Active transcription of rRNA operons condenses the nucleoid in Escherichia coli: Examining the e ect of transcription on nucleoid structure in the absence of transertion[END_REF], Bakshi et al., 2012]. From a point of view where macromolecular crowding is critical for nucleoid compaction, the halt of mRNA synthesis depletes polysomes, and thus results in a lower osmotic pressure on the nucleoid, allowing its expansion. We tested this scenario with our model by using the non-equilibrium steady state shown in Fig. 2.4 A as the initial condition for Eqs. (2.1) to (2.3), switching o mRNA synthesis (α = 0), and integrating forward in time. As shown in Fig. 2.9, the nucleoid expands and spreads over most of the intracellular space. The nucleoid does not take over the entire cell because there are pockets of free ribosomes at both cell poles, which prevent the DNA from occupying these spaces.

The nucleoid reaches its expanded steady state in ∼ 30 min, which is in good agreement with experimental data [START_REF] Cabrera | Active transcription of rRNA operons condenses the nucleoid in Escherichia coli: Examining the e ect of transcription on nucleoid structure in the absence of transertion[END_REF], see Fig. 2.9 B, despite the fact that parameters like total cell length or growth rate used in our model may not match the (unreported) experimental ones. The agreement is due to the fact that the expansion process is primarily driven by the degradation of mRNAs, hence β is the most relevant parameter for this process. In Fig. 2.9 B the bulk of the expansion happens in the rst 10 min -a timescale consistent with the half-life of mRNA (5 min) and with the experimental data. Therefore, the most relevant parameter for the expansion process is the mRNA degradation rate and, given the quantitatively good behaviour of the model, we conjecture that in the experiments by Cabrera and colleagues this parameter takes ). The turquoise data points were obtained from Cabrera et al. [START_REF] Cabrera | Active transcription of rRNA operons condenses the nucleoid in Escherichia coli: Examining the e ect of transcription on nucleoid structure in the absence of transertion[END_REF], for the case of a cell treated with rifampicin (which blocks mRNA transcription).

a similar value to the one used here (which was obtained from Ref. [START_REF] Bernstein | Global analysis of Escherichia coli RNA degradosome function using DNA microarrays[END_REF]).

Halt in translation

It has been shown that, after inhibition of translation with antibiotics such as kasugamycin or chloramphenicol, the nucleoid contracts [START_REF] Bakshi | Time-dependent e ects of transcription-and translation-halting drugs on the spatial distributions of the Escherichia coli chromosome and ribosomes[END_REF]. We therefore tested our model to see if it can explain the contraction of the nucleoid after treatment with kasugamycin. Kasugamycin is known to reduce the abundance of polysomes by inhibiting translation initiation and to substantially compact the nucleoid [START_REF] Xiang | Interconnecting solvent quality, transcription, and chromosome folding in Escherichia coli[END_REF]. Hence, to reproduce the e ect of kasugamycin in the model we set the rate of ribosome binding to polysomes to zero (k on = 0) and integrate foward in time from the non-equilibrium steady state shown in Fig. 2.4 A. The results are reported in Fig. 2.10. There is a signi cant contraction of the nucleoid after ∼ 2 min, which is consistent with the timescale on which ribosomes unbind from polysomes k -1 o = 40 s. Nevertheless, the contraction of the nucleoid in the model is smaller than the one reported in Ref. [START_REF] Xiang | Interconnecting solvent quality, transcription, and chromosome folding in Escherichia coli[END_REF], potentially re ecting the e ect of transertion on nucleoid size [START_REF] Bakshi | Time-dependent e ects of transcription-and translation-halting drugs on the spatial distributions of the Escherichia coli chromosome and ribosomes[END_REF], which is not taken into account in our model. However, the fact that the model predicts a contraction of the nucleoid if ribosomes stop binding to mRNAs suggests that the steric interaction of bare mRNAs with the DNA is enough to compact the nucleoid and that polysomes are not required for nucleoid compaction, which helps explain recent observations [START_REF] Xiang | Interconnecting solvent quality, transcription, and chromosome folding in Escherichia coli[END_REF].

It was also noted that treatment with kasugamicin seems to destabilise the positioning mechanism of daughter nucleoids and cause their coalescence [START_REF] Xiang | Interconnecting solvent quality, transcription, and chromosome folding in Escherichia coli[END_REF]. Given that treatment with kasugamicin decreases the width of the nucleoid [START_REF] Bakshi | Time-dependent e ects of transcription-and translation-halting drugs on the spatial distributions of the Escherichia coli chromosome and ribosomes[END_REF], this will increase the di usion of crowders around the nucleoid, implying that our one-dimensional model may not be valid anymore, in accordance with the ndings of the previous section. Therefore, the positioning mechanisms suggested in this chapter may be less e ective after treatment with kasugamicin and could eventually cause the coalescence of two daughter nucleoids. As argued below, further research should include a more comprehensive description of the radial degrees of freedom in the spatial organisation of E. coli to better understand this kind of perturbations and the general mechanisms of self-organisation.

Discussion of the results

In this chapter, we investigated the physical origins of the intracellular localisation of DNA, messenger RNAs (mRNAs), and ribosomes in bacteria. This is a topic of general interest due to its far-reaching consequences, such as, the spatial organisation of transcription and translation [START_REF] Gray | Nucleoid size scaling and intracellular organization of translation across bacteria[END_REF], Weng et al., 2019], chromosome positioning and segregation [Wu et al., 2019, Joshi et al., 2011], and a wide range of cellular processes regulated by the nucleoid that excludes many macromolecules from the volume which it occupies [START_REF] Coquel | Localization of protein aggregation in Escherichia coli is governed by di usion and nucleoid macromolecular crowding e ect[END_REF].

Based on steric interactions among DNA, mRNAs and ribosomes, we developed a holistic model for the spatial organisation of the bacterial nucleoid, which is intuitive and can quantitatively describe the behaviour of the nucleoid and account for experimentally tested perturbations. The model predicts the formation of a phase-separated nucleoid, whose size is in agreement with experimental measurements [Wu et al., 2019] for cells smaller than 10 µm (Fig. 2.5). Beyond this cell length, our model is no longer accurate, for reasons that may include the lack of connectivity among modeled DNA segments, uncertainties in the concentration of crowders, and molecular components which have not been explicitly incorporated in the model, such as nucleoidassociated proteins [START_REF] Dame | Chromosome organization in bacteria: mechanistic insights into genome structure and function[END_REF] or topoisomerases [START_REF] Stuger | DNA supercoiling by gyrase is linked to nucleoid compaction[END_REF] (see Section 2.1).

Moreover, the model highlights the importance of transcription and translation in regulating the size of the nucleoid. On the one hand, it can account for nucleoid expansion as a result of a halt in mRNA synthesis, demonstrating that the progressive degradation of crowders could be the physical cause of the expansion. Indeed, the timescales on which such expansion happens matches the ones observed experimentally [START_REF] Cabrera | Active transcription of rRNA operons condenses the nucleoid in Escherichia coli: Examining the e ect of transcription on nucleoid structure in the absence of transertion[END_REF], and coincides with the 2.5. D timescales of mRNA turnover. On the other hand, halting translation (by precluding ribosomes from binding to polysomes or mRNAs) further compacts the nucleoid because bare mRNAs are still large enough to sterically compact the nucleoid and the number of free ribosomes is increased.

Our results also underline the role of non-equilibrium e ects in the regulation of nucleoid position. The nucleoid is known to localise at midcell [Wu et al., 2019], and we demonstrate that the synthesis of mRNAs and their expulsion from the nucleoid caused by steric e ects can give rise to this positioning -see Fig. 2.8. In fact, a perturbation from the central position of the nucleoid induces an osmotic-pressure di erence between the two cell poles, which pushes the nucleoid back to midcell. The timescale for this centring depends on both the time it takes to establish an osmotic-pressure di erence, which is set by the mRNA turnover time, and the drag experienced by the nucleoid. This drag may be underestimated in our model, because we do not include e ects that could slow down nucleoid centering, e.g. the transient attachment of the nucleoid to the membrane by proteins that are simultaneously being transcribed, translated, and inserted into the membrane, also known as transertion [START_REF] Gorle | DNA condensation in live E. coli provides evidence for transertion[END_REF]. Furthermore, our model shows that non-equilibrium e ects are responsible for the ubiquitous nucleoid splitting and localisation at 1/4 and 3/4 positions along the long cell axis. Indeed, the synthesis of mRNAs within the nucleoid, without additional active processes, is a robust mechanism to make the daughter nucleoids localise at 1/4 and 3/4 positions, as observed experimentally [Wu et al., 2019].

Turing patterns at the cellular scale [START_REF] Wu | Multistability and dynamic transitions of intracellular Min protein patterns[END_REF]] also display features, such as, the emergence of a characteristic lengthscale, that could seem similar to the ones produced by our model, in view of Fig. 2.11. These non-equilibrium patterns have been used to investigate many biological features on a cellular and sub-cellular scale, such as the positioning of protein clusters in E. coli [START_REF] Murray | Self-organization and positioning of bacterial protein clusters[END_REF]. However, unlike Turing patterns, our model predicts a phase separation due exclusively to steric interactions, and in the absence of non-equilibrium e ects, see Fig. 2.7 A, upper pannel. The patterns reported here are more closely related to the ones produced by other models with an interplay between an equilibrium free energy and non-equilibrium dynamics, e.g., the formation of FtsZ rings in bacteria [START_REF] Shlomovitz | Membrane-mediated interactions drive the condensation and coalescence of FtsZ rings[END_REF], or other nonequilibrium phase-separation models [START_REF] Li | Non-equilibrium phase separation with reactions: a canonical model and its behaviour[END_REF], such as models of growing droplets [START_REF] Zwicker | Growth and division of active droplets provides a model for protocells[END_REF]. Nevertheless, our model provides a conceptually simpler framework to produce these patterns. In fact, unlike a model of physically growing droplets, our analysis involves a conserved order parameter -the total number of DNA segments -and the e ect of non-equilibrium terms -mRNA production and degradation -is limited to nucleoid reshaping and repositioning. Despite its simplicity, our model produces a number of experimentally observed patterning e ects, such as nucleoid centering at midcell, splitting and positioning of sister lobes during cell division.

In addition, this patterning is not limited to nucleoid splitting into two sister lobes, because our model predicts that the nucleoid can split into more than two lobes, whose size is given by a characteristic length and whose position is tightly controlled, as can be seen in Fig. 2.11. Experimentally, in long lamentously growing cells (where cell division is inhibited but DNA replication is allowed to continue), nucleoids are observed at tightly controlled positions and distances [START_REF] Wehrens | Size laws and division ring dynamics in lamentous Escherichia coli cells[END_REF], forming a very similar pattern to the one obtained here, albeit with a shorter characteristic length. In Ref. [START_REF] Wehrens | Size laws and division ring dynamics in lamentous Escherichia coli cells[END_REF] a separate nucleoid appears every To enable comparison between di erent panels, the lengthscale in all panels was kept constant. As in Fig. 2.7, for cells larger than ∼ 8 µm the nucleoid splits in two lobes. Furthermore, for cells around ∼ 16 µm the nucleoid has three distinct lobes, suggesting that this may be a pattern with a characteristic length that exists also for longer cells.

∼ 2.25 µm the cell grows in length. This value is roughly half of the one predicted by this model, possibly due to uncertainties in the parameters such as the di usion coe cient of mRNAs within the nucleoid.

Given its generality, our analysis is not necessarily restricted to the nucleoid of prokaryotic cells [START_REF] Cohan | Making the case for disordered proteins and biomolecular condensates in bacteria[END_REF]. For instance, division of certain phase-separated condensates has been experimentally related to non-equilibrium processes, as is the case for the ParABS partition system, which creates phase-separated condensates of DNA and ParB around parS sites, whose division is controlled by the activity of ParB's ATPase activity on ParA [START_REF] Guilhas | ATP-driven separation of liquid phase condensates in bacteria[END_REF]. The activity-driven nucleoid division described in our model may thus constitute a general strategy employed by cells to control the structure and positioning of membraneless compartments.

Future directions: Extensions and experimental tests

Our work implies that steric interactions make the bacterial cytoplasm an e ectively poor solvent for the chromosome, as recently indicated by experiments [START_REF] Xiang | Interconnecting solvent quality, transcription, and chromosome folding in Escherichia coli[END_REF]. However, steric interactions may not be the only contribution to the poor-solvent quality of the cytoplasm. Other types of intermolecular interactions [Odijk, 1998] or the e ect of nucleoid-associated proteins [START_REF] Dame | Chromosome organization in bacteria: mechanistic insights into genome structure and function[END_REF] could also a ect the solvent quality of the cytoplasm and, therefore, the organisation of the nucleoid in the cell. Notably, certain proteins or enzymes, such as Dps [START_REF] Janissen | Global DNA compaction in stationary-phase bacteria does not a ect transcription[END_REF] and DNA gyrase [START_REF] Stuger | DNA supercoiling by gyrase is linked to nucleoid compaction[END_REF], have been shown to modify the compaction of the nucleoid. Given that the radii of the steric interaction might be an e ective one, phenomenologically, one could account for perturbations where these molecules are involved by varying the parameter ρ (see Section 2.2.1) and reduce it to account for their presence. Still, for future studies, both theoretical and experimental, research into these other regulators of the nucleoid size could yield a more complete picture of its organisation, and improve the accuracy of the results presented here.

Another of the limitations of this study is the reduction from the three dimensions of an E. coli cell to the single one considered here. In reducing the system to one dimension, we tacitly assumed that the nucleoid occupies the whole cross-section. While normally the nucleoid does spread over most of the cross-section, upon treatment by kasugamycin and chloramphenicol the width of the nucleoid decreases [START_REF] Bakshi | Time-dependent e ects of transcription-and translation-halting drugs on the spatial distributions of the Escherichia coli chromosome and ribosomes[END_REF] with two important consequences: Since the nucleoid does not occupy the whole cross-section, the osmotic pressures in di erent regions of the cells can be equilibrated without the need for lateral displacement of the nucleoid, which removes the mechanism that positions the nucleoid at midcell (or 1/4 and 3/4 for daughter nucleoids). In addition, if the nucleoid does not occupy the whole cross-section of the cell, mRNAs might be able to escape radially from the nucleoid, which can induce a deviation from the scaling found in Eq. (2.49) for the length at which the nucleoid splits. Indeed, previously separated daughter nucleoids have been observed to coalesce upon treatment with kasugamycin and chloramphenicol [Schaechter andLaing, 1961, Xiang et al., 2021] and we hypothesise that the decrease in nucleoid width induced by these drugs destabilises the positioning of daughter nucleoids. Then, if the nucleoids come into contact due to uctuations, they are driven to coalesce by short-range depletion interactions due to ribosomes, mRNAs, and any remaining polysomes. Therefore, we expect that exploring the three-dimensional dynamics of the system might explain C 2. O several phenomena that for the moment remain unaccounted for.

Moreover, the present work has been devoted to the physical principles governing the global organisation and positioning of the bacterial nucleoid. However, the nucleoid also has a rich internal dynamics [START_REF] Fisher | Four-dimensional imaging of E. coli nucleoid organization and dynamics in living cells[END_REF]. Given that the DNA plectoneme persistence length is close to the cell radius, future studies may investigate whether orientational order emerges, and whether it has implications for intra-nucleoid organisation.

As outlined above, our model (like every model) has certain limitations. Therefore, an important part of building a theory is knowing when it will collapse and will be no longer accurate or able to explain the observations. The grounds for this testing must also be laid.

In the case of the present theory for the organisation of the bacterial nucleoid, the onedimensional assumption, which turns out to be crucial for nucleoid positioning, is one its limits, as mentioned above. Yet, its consequences are not restricted to the coalescence of daughter nucleoids, as a substantial decrease in the radius of the nucleoid should enable the di usion of macromolecules around it, destabilising too the centring mechanism. Therefore, if the nucleoid shrinks in the radial dimension, while the average position of the nucleoid may still be the centre4 , there should be larger uctuations in the position of the nucleoid. These should be measurable by the second moment of the distribution of the nucleoid position, if measurements of enough resolution can be made for a large enough number of cells and long enough times, and treatment by kasugamycin (or erythromycin [START_REF] Sanamrad | Single-particle tracking reveals that free ribosomal subunits are not excluded from the Escherichia coli nucleoid[END_REF]) should result on an increase of the variance of the position.

However, treating the cells with kasugamycin induces a large perturbation in the metabolism of a cell, since not only does it make the nucleoid shrink but it also halts translation, with the widespread consequences this may carry. Thus, alternative approaches that may yield similar outcomes should also be envisaged. In particular, quantitative measurements of the nucleoid in stationary phase E. coli could be a good candidate.

When nutrients are scarce and cannot sustain exponential growth, E. coli cells enter stationary phase and stop reproducing [Nyström, 2004, Santos et al., 2005]. Upon cessation of growth, the synthesis of Dps is induced [START_REF] Almiron | A novel DNA-binding protein with regulatory and protective roles in starved Escherichia coli[END_REF], which compacts the nucleoid at least in the longitudinal direction [START_REF] Janissen | Global DNA compaction in stationary-phase bacteria does not a ect transcription[END_REF]. If a radial contraction of the nucleoid was also measured, then, according to our theory, there should be a destabilising e ect in the central location of the nucleoid in stationary phase (larger dispersion of the measured positions even if the average is still at midcell). Moreover, stationary phase cells are often associated with a slower transcriptomic activity, which should further amplify the instability of nucleoid location.

Therefore, quantitative assessment of the stationary-phase E. coli nucleoid compaction and positioning could also be used to compare with the theory exposed here and explore its limitations. However, care should be taken when extrapolating the results obtained here to the stationaryphase nucleoid, since the highly expressed Dps protein would increase the attractive interactions between DNA segments (in a phenomenon akin to other protein-DNA co-condensations [START_REF] Quail | Force generation by protein-DNA co-condensation[END_REF]) and, eventually, it may outcompete repulsive interactions, thus invalidating the picture based on excluded-volume interactions proposed here [START_REF] Andersen | Relationship between the hard-sphere uid and uids with realistic repulsive forces[END_REF] The examples exposed in the rst part of this thesis highlight the importance of interactions between solutes in the cell cytoplasm. These interactions give rise to emergent behaviour like phase-separated organelles (e.g. the nucleoid in certain bacterial cells, see Chapter 2) or control the size of protein aggregates, as found in the rst chapter. Nonetheless, Part I was devoted to two particular examples and one would like to have a theory as general as possible for a wide variety of processes, to increase our understanding of them and be able to distinguish between the generic properties of these mechanisms and the details of each one of them.

An umbrella under which many of these processes can be described is that of reaction-di usion equations. More precisely, we will consider a network of chemical reactions that is driven out of equilibrium to represent the reactions that may happen in the cytoplasm of a cell. In addition, di usion must be taken into account if we aim at explaining the spatial organisation within a cell. However, unlike classical reaction-di usion systems, which are assumed to be ideal (i.e. solutions with non-interacting solutes), we need to consider interactions between the components of the solution to better account for the complex dynamics of the crowded cytoplasm and, thus, build more realistic models that describe better the physical reality of a cell. Altogether, the framework we aim at building may be termed non-ideal reaction-di usion systems.

Beyond the examples of Part I, it is becoming increasingly clear that, in general, the cytoplasm of a cell does not behave like an ideal solution [START_REF] Zielinski | Nonideal solute chemical potential equation and the validity of the grouped solute approach for intracellular solution thermodynamics[END_REF]. There is a plethora of interactions between the solutes that take place within the cytoplasm, mostly involving proteins and other macromolecules, but also ions. Some of the most common interactions that are relevant in the cellular cytoplasm are steric or crowding e ects -see Part I and Refs. [START_REF] Zhou | Macromolecular crowding and con nement: Biochemical, biophysical, and potential physiological consequences[END_REF], Mitchison, 2019] -and electrostatic interactions (e.g. [Fall andKeizer, 2001, Wang et al., 2018]). Arguably, the most striking phenomenon caused by these interactions is the emergence of phaseseparated condensates within the cytoplasm, also known as membraneless organelles particularly in the cell biology literature, which are now widely studied [START_REF] Li | Phase transitions in the assembly of multivalent signalling proteins[END_REF], Su et al., 2016, Brangwynne et al., 2009]. The composition of these membraneless organelles is di erent to that of the cytoplasm, as they are typically enriched in a certain type of molecules while they exclude others [START_REF] Ditlev | Who's in and who's out-compositional control of biomolecular condensates[END_REF]. Moreover, it has been hypothesised that they spatially control biochemical reactions [START_REF] Banani | Biomolecular condensates: organizers of cellular biochemistry[END_REF], by modulating their rates and speci city within the condensate. In this sense, the nucleoid described in Chapter 2 shares all of the properties described here for membraneless organelles, since it controls the spatial distribution of transcription and translation, and it excludes large macromolecules from its volume.

Therefore, phase separation phenomena play an important role in the internal spatial organisation of cells and their regulation is crucial for many cellular functions. One of the ways cells can dynamically control the appearance of membraneless organelles (and its composition or function) is through chemical reactions, notably post-translational modi cations like phosphorylation [START_REF] Tsang | Phosphoregulated FMRP phase separation models activity-dependent translation through bidirectional control of mRNA granule formation[END_REF], Kim et al., 2019] or methylation [START_REF] Nott | Phase transition of a disordered nuage protein generates environmentally responsive membraneless organelles[END_REF]. However, phase separation is also triggered by changes in the environment [START_REF] Franzmann | Phase separation of a yeast prion protein promotes cellular tness[END_REF], Lyon et al., 2021], establishing them as potential sensing and regulatory mechanisms.

From the theoretical perspective, the interplay between interactions within the solution and non-equilibrium chemical reactions has also been widely studied. However, most of these e orts [Huberman, 1976, Glotzer et al., 1995, Wurtz and Lee, 2018, Li and Cates, 2020] were based on e ective reaction-di usion models that, while they can describe patterning and other nonequilibrium phenomena in a simple way, they lack thermodynamic consistency, as they treat chemical reactions as ideal -by modelling their dynamics with mass-action kinetics (MAK)while the interaction-in uenced di usion that drives phase-separation is not, and there is no explicit non-equilibrium driving. Conversely, a thermodynamically consistent description would require that, in the same way in which di usion is governed by a free energy (that takes into account the interactions), the dynamics of chemical reactions must also re ect this free-energetic dependency. Early progress at addressing non-ideal reaction-di usion systems in a thermodynamically consistent model was limited to a linear stability analysis for the dynamics of binary solutions [START_REF] Carati | Chemical freezing of phase separation in immiscible binary mixtures[END_REF]]. More recently, some works have aimed at establishing a deterministic theory for non-ideal chemical reaction networks (CRNs) [START_REF] Avanzini | Nonequilibrium thermodynamics of non-ideal chemical reaction networks[END_REF] and exploring minimal examples for pattern formation [Bazant, 2013, Kirschbaum andZwicker, 2021] but the link between non-equilibrium CRNs and phase-separation has not yet been elucidated in full generality.

In this part, we aim at building a thermodynamically consistent framework for interacting reaction-di usion systems that, thus, may exhibit phase-separation at steady state. Therefore, a complete theory of non-ideal CRNs is necessary and here previous e orts are complemented by analysing the behaviour of non-ideal CRNs in the stochastic limit and exploring the consequences of the topology of the network. We do so by rst constructing a framework which satis es that, in the absence of explicit non-equilibrium driving, the system relaxes to thermodynamical equilibrium (Chapter 3). This enables us to naturally adapt and generalise results from the well-established theory of ideal CRNs, explore the connection between non-equilibrium CRNs and phase separation, and obtain Lyapunov functionals for a particular class of CRNs, known as complex-balanced (Chapter 4).

A. M. Turing already stated the di culty of building a general theory for patterning processes and, here, we will fail in doing so, as the applicability of most of our results is limited to complex-balanced networks. However, important consequences can be derived from these types of networks, see Chapter 5, where we also discuss on the implications of considering more general networks.

Chapter 3

Introduction to the description and thermodynamics of Chemical Reaction Networks

We will now take the approach of the formal kineticist 1 and de ne the chemically reacting systems in terms of elementary reactions. An elementary reaction is a chemical reaction with no hidden intermediate steps [Kondepundi and Prigogine, 2015], that is, the most basic unit of a chemical transformation 2 , which is de ned by a set of stoichiometric coe cients and a rule relating reaction rates to composition (as well as other physical variables), which for ideal solutions usually is the law of mass action kinetics [START_REF] Horn | General mass action kinetics[END_REF]. But, most importantly, the driving force of elementary reactions (whether in or out of thermodynamic equilibrium) is their free energies of reaction or a nity (no matter how complex these free energies may be), provided a sort of local equilibrium can be de ned [Kondepundi and Prigogine, 2015].

Therefore, the thermodynamically consistent approach to chemical reaction networks (CRNs) in non-ideal solutions must be built upon the concept of elementary reactions because it is only by taking into account every process happening in the system that the laws of themodynamics can be recovered for CRNs, even when these are driven out of thermodynamic equilibrium. Once the dynamics of a CRN have been written down in a thermodynamically consistent way, we will be able to properly de ne thermodynamic quantities, such as entropy production or work; and extend to non-ideal solutions results that are well-known for ideal CRNs.

In order to construct a thermodynamically consistent framework for non-ideal CRNs we will impose two conditions:

• A CRN must relax to thermodynamic equilibrium in the absence of external work.

1 "The formal kineticist, on the other hand, takes a macroscopic viewpoint and his primitive concept is the elementary reaction." F. Horn and R. Jackson, General Mass Action Kinetics, 1972. 2 Accordig to the IUPAC's Gold Book, an elementary reaction is "a reaction for which no reaction intermediates have been detected or need to be postulated in order to describe the chemical reaction on a molecular scale. An elementary reaction is assumed to occur in a single step and to pass through a single transition state." In practice, however, the timescales of relaxation of any additional transition states and any intermediate metastable states must be taken into account for a practical but meaningful de nition of elementary reaction.

61 62 C

D CRN

• The external work done on the CRN must be speci ed explicitly.

Here we will restrict ourselves to CRNs where a subset of chemical species is chemostatted, that is, kept at a certain chemical potential as a consequence of external work. However, more complex procedures can be devised to keep a CRN out of thermodynamic equilibrium.

Finally, in the rest of this work, we will use the word equilibrium to refer to thermodynamic equilibrium. Conversely, we will call steady state the state of a system with vanishing time derivatives3 .

Description of Chemical Reaction Networks

We start by describing a CRN comprising N di erent chemical species and M di erent reversible reaction pathways. The requirement that every reaction is reversible is a requirement for the thermodynamic consistency of our description of the CRN if all of the reaction a nities are nite. In full generality, a reaction ρ belonging to the CRN can be speci ed as follows:

i r ρ i X i i s ρ i X i , (3.1)
where X i denotes the di erent species labeled with the subindex i. In the rest of the document the indices i and j will be reserved for chemical species. The matrix r ρ i denotes the stoichiometry of the forward reaction (speci es number of reactants of type i required by the reaction ρ) and s ρ i that of the backwards reaction (products of type i produced by the reaction ρ). Note that, given that the reactions are taken to be reversible, the distinction between reactants and products is arbitrary. Therefore, v ρ i = s ρ i -r ρ i is the net amount of i species created or destroyed along the forward direction of the reaction ρ (denoted as +ρ). We also de ne the vectors

v ρ = (v ρ 1 , • • • , v ρ N ), r ρ = (r ρ 1 , • • • , r ρ N ) and s ρ = (s ρ 1 , • • • , s ρ N )
for compactness, which are vectors carrying the information of the type and number of species involved in ρ as net change, reactants or products, respectively. Even more compact, we de ne the matrices

V = (v 1 , • • • , v M ), R = (r 1 , • • • , r M ) and S = (s 1 , • • • , s M )
, whose elements are v ρ α , r ρ α and s ρ α , respectively. We nally de ne the notion of complex z, which is the set (number and type) of particles that take part in a chemical reaction as reactants (z = r ρ ) or products (z = s ρ ). A single complex z may appear in more than one reaction within the network.

Example 1 Let us consider the following CRN:

A + B C, B D, (3.2)
whose stoichiometric matrices for the forward and backward reactions are

R = 1 1 0 0 0 1 0 0 , S = 0 0 1 0 0 0 0 1 , (3.3) 3.1. D CRN 63
where each column denotes a particular species and each row a particular reaction of our CRN. Each of the rows in these two matrices corresponds to the vectors r ρ (for the matrix R) and S ρ (for the matrix s).

There are, therefore, four complexes in this network, two reactant complexes and two product complexes.

Stochastic Chemical Master Equation

If the solutes in a dilute4 solution di use fast (with respect to the typical timescale of chemical reactions) the system may be viewed as well-mixed and can, therefore, be described in terms of a single homogeneous concentration of each of the species across the entire system. In this case a state of the system (number and type of particles) is completely determined by the vector

n = (n 1 , • • • , n N )
where n i is the number of particles of type i. CRNs where the copy numbers of certain species are low can show large uctuations and high stochasticity. To describe a stochastic system, each of the states of the system n will have an associated probability measure P (n, t) at each instant of time t. The dynamics for the probabilities of states of homogeneous CRNs is given by the Chemical Master Equation (CME) which we write explicitly taking into account the reversibility of all reactions [Gillespie, 1992]:

∂P (n, t) ∂t = ρ f +ρ (n -v ρ )P (n -v ρ ) + ρ f -ρ (n + v ρ )P (n + v ρ ) - ρ [f +ρ (n) + f -ρ (n)]P (n),
(3.4) where the summations over ρ run over all the reactions in the CRN and the rate of the transitions in the network is given by the propensity function f ρ . In ideal dilute solutions and for elementary reactions, the propensity function f ρ takes the form

f +ρ (n) = k +ρ i n i ! (n i -r ρ i )!
, and

f -ρ (n) = k -ρ i n i ! (n i -s ρ i )! , (3.5)
which is known as mass-action kinetics (MAK) for the stochastic CME. In the following section we develop a generalisation of this propensity function for non-ideal systems based on the detailedbalance condition.

In most cases, the dynamics of the CRN will tend to a stationary distribution in the long-time limit. If this is the case, the steady state of the dynamics of the network will be de ned by a zero time derivative of P (n, t) and will therefore satisfy the following relation:

ρ f +ρ (n -r ρ )P (n -r ρ ) + ρ f -ρ (n -s ρ )P (n -s ρ ) = ρ (f +ρ (n) + f -ρ (n))P (n). (3.6)
In general, it is not possible to obtain an analytical form for the steady-state distribution of the CME. However, for a broad class of steady states its probability distribution can be constructed analytically, as it will be shown in the following chapter.

C

3. D CRN

Deterministic description of a Chemical Reaction Network

For large particle numbers, by taking the averages of eq. ( 3.4) and making mean-eld assumptions -that is, assuming vanishing correlations -one can derive a set of equations for the concentrations in the macroscopic limit, that is, when n i /V = c i , for both n i and V large. T. G. Kurtz [Kurtz, 1972] showed that these two approaches are indeed equivalent in the large volume limit, at least for MAK. In this limit, the state of the system is fully speci ed by a set of concentrations (c 1 , • • • , c N ) and one obtains the following classical set of equations for the dynamics of the concentrations c i in a CRN [START_REF] Schnoerr | Approximation and inference methods for stochastic biochemical kinetics-a tutorial review[END_REF]:

∂c i ∂t = ρ v ρ i (J +ρ -J -ρ ), (3.7)
where the currents J still need to be determined. For the classical case of MAK, in the deterministic description, the currents take the form

J +ρ = k +ρ i c r ρ i i , and J -ρ = k -ρ i c s ρ i i . (3.8)
While both the deterministic and stochastic descriptions refer to the same system, the former one is only accurate for large particle numbers, also known as the thermodynamic limit, where uctuations are negligible.

Thermodynamical constraints on the dynamics of Chemical Reaction Networks

In the previous section we introduced the general description of CRNs, both at the stochastic and deterministic level, but, in both cases, a choice for the propensity functions or currents must be made. For ideal solutions, the most common choice is MAK, as outlined above. However, here we are considering solutes that interact among themselves, which are, therefore, not ideal. In this section, our aim is to specify the propensity functions or currents for non-ideal systems by assuming that, if there is no explicit non-equilibrium driving, the system must relax to thermodynamic equilibrium.

Stochastic dynamics

At thermodynamic equilibrium, detailed balance must necessarily hold for every reaction ρ, in order to ful ll the second law of thermodynamics. Detailed balance implies that the probability ux across a reaction ρ in the forward direction must equal the probability ux in the backward direction, thus preserving time reversal symmetry in the system and precluding from systematically extracting work from the uctuations in the system [Sekimoto, 2010]. Mathematically, for stochastic systems, detailed balance takes the following form

P eq (n)f +ρ (n) = f -ρ (n + v ρ )P eq (n + v ρ ), (3.9) 3.2. T CRN 65
where the equilibrium probability distribution P eq (n) for closed stochastic systems (no energy or matter5 exchanged with the environment) is given by the canonical Boltzmann distribution:

P eq (n) = 1 Z e -βF (n) , (3.10)
β being the inverse of temperature times the Boltzmann constant, F (n) the Helmholtz free energy of the system at state n and Z a normalisation factor, known as the partition function in the statistical physics literature. Importantly, the free energy F that appears in Eq. (3.10) is not necessarily an ideal free energy but can be as complex as needed and include interactions between the di erent constituents of the system.

Detailed balance, as expressed in Eq. (3.9), yields the following constraint for the propensity functions: .11) implying that the reaction rates do depend on the interactions within the system through the changes in free energy provoked by the reaction. Then, we can choose the following functional form for the propensities f ρ :

f +ρ (n) f -ρ (n + v ρ ) = P eq (n + v ρ ) P eq (n) = e -β[F (n+v ρ )-F (n)] , ( 3
f +ρ (n) = k ρ e β[F (n)-F (n-r ρ )] , f -ρ (n + v ρ ) = k ρ e β[F (n+v ρ )-F (n+v ρ -s ρ )] , (3.12)
where k ρ is the reaction constant, which has to be equal for both the forward and the backward reaction6 . This is a di erence with the usual MAK description, where the dependency of the rates on the standard-state chemical potential of the species (the part of the chemical potential that does not depend on the state of the system n) is hidden in the rate constants [Kondepundi andPrigogine, 2015, Rao andEsposito, 2016]. This form for the propensity functions satis es Eq. (3.11) and means the rate of the reactions is a function of the free energy of the reactant complex.

Note that, while we impose detailed-balance conditions to the rates at equilibrium, this does not mean they are only valid for the equilibrium state, it only ensures that the system relaxes to equilibrium in the absence of external driving or work. By doing so, we found a functional form for the propensity functions (3.12) that we will consider to be valid even if the system is not at thermodynamic equilibrium. This is inspired by the fact that rates of the form of MAK are thought to be an accurate description of reactions both at equilibrium and far away from equilibrium [START_REF] Groot | Non-equilibrium Thermodynamics[END_REF].

Mass Action Kinetics can be derived from an ideal free energy

The form (3.12) for the propensity functions is not unique but is particularly appealing given that it reduces to MAK for ideal systems. In an ideal system, at least in a lattice model for the solution C

3. D CRN (see Appendix C.1) with total volume V = i n i , the free energy is

F id = i n i µ 0 i + β -1 i log(n i !) -log(V !) , (3.13)
where µ 0 i is the standard-state chemical potential of species i. Then, the rates take the following form:

f +ρ (n) = k ρ e β[F (n)-F (n-r ρ )] = k ρ e β i r ρ i µ 0 i (V -i r ρ i )! V ! i n i ! (n i -r ρ i )! , (3.14)
where (Vi r ρ i )!/V ! can be approximated by V -i r ρ i . After rede ning

k +ρ = k ρ exp(β i r ρ i µ 0 i ) V i r ρ i , (3.15)
we nd that the propensity functions (3.12) for an ideal free energy matches the propensity function for stochastic MAK, Eq. (3.5). Note that, after rede ning the reaction constant k +ρ , the forward and backward reaction constants are no longer equal.

Example 2 Let us consider the second reaction in Eq. ( 3.2):

B D (3.16)
Then the propensity functions, according to Eq. ( 3.12), are

f +ρ (n) = k ρ exp(β[F (n) -F (n -r ρ )]) f -ρ (n) = k ρ exp(β[F (n) -F (n -s ρ )]).
(3.17)

Let us now de ne the chemical potential of the i-th species for stochastic (small) systems as follows:

µ i (n) = F (n) -F (n -e i ), (3.18)
where e i is the unit vector with all entries zero except for the i-th entry that is 1. Then, the propensity functions of unimolecular reactions (like B D) can be written as

f +ρ (n) = k ρ exp[βµ B (n)] f -ρ (n) = k ρ exp[βµ D (n)], (3.19)
which, again, reduces to mass action kinetics for unimolecular reactions if we consider an ideal system.

Deterministic dynamics

In systems with large particle numbers, where a deterministic description, as the one presented in Section 3.1.2, is accurate, one can adapt the rates (3.12) obtained for the stochastic system to this macroscopic limit. When the particle numbers n are large the di erences in free energies that appear in the rates (3.12) can be expanded in Taylor series and rewritten as (3.20) where µ i is the standard chemical potential of species i, µ i = ∂f (c)/∂c i , and f (c) is the free energy of the system per unit volume in the deterministic notation. Eq. (3.20) becomes exact in the limit when n → ∞, while r ρ remains nite.

F (n) -F (n -r ρ ) i r ρ i µ i ,

T CRN 67

Therefore, the currents in a deterministic and non-ideal CRN, whose dynamics are given by Eq. (3.7), can be written as

J +ρ = k ρ e β i r ρ i µ i , J -ρ = k ρ e β i s ρ i µ i , (3.21)
which is an expression conceptually similar to those given by other approaches to build thermodynamically consistent descriptions for deterministic non-ideal CRNs [START_REF] Avanzini | Nonequilibrium thermodynamics of non-ideal chemical reaction networks[END_REF], Bazant, 2013].

Once again, the currents (3.21) match their ideal MAK counterpart (3.8) if the chemical potentials used in the rates are those of an ideal solution, i.e., µ i = β -1 log c i + µ 0 i . In this case, we obtain, e.g. for the forward reaction, (3.22) which matches the deterministic MAK, Eq. (3.8), and where k +ρ = k ρ exp(β i r ρ i µ 0 i ). In Sections 3.2.1 and 3.2.2, we have written the chemical reaction rates for both stochastic and deterministic systems in terms of their free energies. Given that one of the aims of the work presented here is to describe CRNs in phase-separating systems, in Appendix C.1 we brie y review the Flory-Huggings or regular solution theory with interactions, widely used to model phase separation, in order to obtain free energies and chemical reaction rates (stochastic and deterministic) for phase-separating solutions, results that will be used in the following chapter. In Appendix C.2, we show how to obtain the non-ideal chemical reaction rates from activity coe cients, more common in the chemistry literature.

J +ρ = k ρ e β i r ρ i µ i = k +ρ i c r ρ i i ,

On the rate constant

Until now we have assumed k ρ to be a constant of the reaction and completely insensitive to any change of state in the system. However, in full generality, this constant could depend on the state of the system, as the system is part of the environment the reaction takes place in. These e ects can be disregarded for most cases in ideal or dilute solutions, but in non-ideal ones these e ects may not be negligible, particularly in the case of phase separation, where there could be two (or more) very di erent environments where the chemical reaction takes place, thereby accelerating it or slowing it down. Nevertheless, detailed balance, Eq. (3.11), must continue to hold. This means that the forward reaction constant for a state n must be equal to the backward reaction constant for a state n + v ρ . One way to ensure this equality while keeping the state-dependency of the reaction constants is to make k ρ a function of the state after subtracting the reactant complex, that is nr ρ for the forward case and n + v ρs ρ for the backward one. Noting that the detailedbalance equality is still satis ed by virtue of nr ρ = n + v ρs ρ , we can write down the following more general rates:

f +ρ (n) = k ρ (n -r ρ )e β(F (n)-F (n-r ρ )) , f -ρ (n + v ρ ) = k ρ (n + v ρ -s ρ )e β(F (n+v ρ )-F (n+v ρ -s ρ )) .
(3.23) From a physical perspective, this means that the rate constant can only be dynamically a ected by the entire system except for the complex that takes part in the reaction. The horizontal dimension is the reaction coordinate and the vertical one speci es the height of the reaction free energy landscape F 0 (which is equal to F except for the fact that it does not include the entropic term i log n i !, as is customary [Kondepundi and Prigogine, 2015]). For each phase there are two minima in the free energy landscape corresponding to whether the reaction has occurred (right-hand minimum) or not (left-hand one). The free energy of the system (F ) speci es the free energy of the system before the reaction happens [F 0 (n)] and after [F 0 (n + v ρ )] but it does not determine the height of the reaction free energy barrier ∆F * . Thus, how ∆F * depends on the environement can be crucial for the behaviour of the CRN. The subscripts specify to which of the phases the value of the variable refers to.

In analogy with the classical transition state theory, we can think of the microscopic mechanism of a reaction as a di usion process in a reaction free-energy landscape in the reaction coordinates [Kondepundi andPrigogine, 2015, Hänggi et al., 1990], see Fig. 3.1. Then, the value of the rate constant depends on the height of the free energy barrier ∆F * of the reaction. While the free energies of the reactants and products (the stable states in the reaction landscape) have free energies de ned by F , this is no longer the case for the height of the barrier ∆F * . For systems where more than one phase coexist, how the height of the barrier (and, thus, the value of k ρ ) depends on which phase the reaction takes place in can be crucial for the behaviour of the CRN. Fig. 3.1 graphically summarises the arguments given in this section and it highlights the e ect a di erent phase or environment can have in the reaction free-energy landscapes and, thus, on the reaction constants k ρ . Note that, setting a value of the barrier imposes a relation between the forward and backward reactions and is therefore equivalent to requiring detailed balance to hold at equilibrium.

Rates for non-equilibrium systems

Until now we have considered the rates or propensity functions of chemical reactions that, even if they are not at thermodynamic equilibrium, they will eventually relax to such a state as a consequence of the functional form of these rates. However, we are not interested in systems that relax to equilibrium but in systems that, like living beings, dissipate free-energy to maintain themselves away from thermodynamic equilibrium.

Let us consider that a certain subset of the N species of the system are connected to di erent particle reservoirs (chemostats) that will keep their concentration constant. Then, in general, the 3.3. T C R N 69 system will not relax to equilibrium due to the work done on the system by the chemostats. In this case, the space of possible states of the system is reduced to N dimensions (as the chemostatted species have their concentrations xed) and the thermodynamic constraint on the propensity functions, Eq. (3.9), takes the following form:

f +ρ (n) f -ρ (n + v ρ ) = e -β(F (n+v ρ )-F (n)+ j v ρ j µ j ) (3.24)
where, now, n is a vector that only encodes the particle number of the N non-chemostatted species, F in the rst term in the exponential also refers to the free energy of non-chemostatted species and the summation over j runs only over the chemostatted species. The second term ( j v ρ j µ j ) speci es the chemical work done by the chemostats (with chemical potentials xed at µ j ) when a reaction ρ occurs, which pushes the system out of equilibrium. Here, an approximation has been made by which the chemostatted species are ideal (no interactions with the rest of the chemical species) and abundant, hence, their free energy change can be replaced by the sum of chemical potentials, as is done for deterministic CRNs, see Eq. (3.20), and we are able to neglect uctuations around their chemostatted values. Therefore, the system is now fully speci ed by a vector n that includes the particle numbers of the non-chemostatted species only. Finally, as discussed in Section 3.2.3, the rate constants may have state dependencies and the rates including the contributions from the chemostats can take the following form

f +ρ (n) = kρ g ρ (n -r ρ )e β[F (n)-F (n-r ρ )+ j r ρ j µ j ] , f -ρ (n + v ρ ) = kρ g ρ (n + v ρ -s ρ )e β[F (n+v ρ )-F (n+v ρ -s ρ )+ j s ρ j µ j ] ,
(3.25)

where j r ρ j µ j and j s ρ j µ j include the contribution of the chemostats for the forward and backward reactions, respectively. For the sake of clarity, in the remainder of the text, we will reserve the index j to denote the chemostatted species and the index i to denote the rest of chemical species (non-chemostatted). Finally, kρ is constant and any state dependency k ρ may have is absorbed into the function g ρ . We note that a similar, and equally general, form for the kinetics has been used in other contexts, such as the modelling of molecular motors [START_REF] Jülicher | Modeling molecular motors[END_REF]. Analogously, in the deterministic limit, the above propensity functions become the following currents

J +ρ (c) = kρ g ρ (c)e β[ i r ρ i µ i + j r ρ j µ j ] , J -ρ (c) = kρ g ρ (c)e β[ i s ρ i µ i + j s ρ j µ j ] , (3.26)
where, µ i is the chemical potential of the non-chemostatted species and is thus a dynamic variable, while µ j corresponds to the chemical potential of chemostatted species and is therefore xed.

Thermodynamics of Chemical Reaction Networks

Until now we have mostly considered the dynamics of non-ideal CRNs, but for the rest of the chapter we will focus on their energetics and review certain thermodynamic relationships valid for stochastic systems far away from equilibrium. These relationships will be useful in the following chapter, as they will shed light into the thermodynamic interpretations of the results obtained there.

C 3. D CRN

Thermodynamics of closed Chemical Reaction Networks

In this section we consider a closed stochastic non-ideal CRN. This implies that the system, if prepared in an arbitrary state, will always relax to its equilibrium state, since it is not connected to any chemostat or any other source of work. Therefore, we can de ne a thermodynamic potential [START_REF] Rao | Conservation laws and work uctuation relations in chemical reaction networks[END_REF]]

G(n, t) = k B T log P (n, t) + F (n) (3.27)
whose average value, taken over the probabilities P (n, t), acts as a Lyapunov function of the system and takes its minimum value at thermodynamical equilibrium. Thus, G can be called the stochastic free energy of the system, as it includes a contribution from the probability P of the system being at state n and time t.

Then, entropy of the system is [START_REF] Schmiedl | Stochastic thermodynamics of chemical reaction networks[END_REF]]

S(n, t) = - ∂G(n, t) ∂T = -k B log P (n, t) + s(n), (3.28)
where

s(n) = -∂F (n)/∂T . The enthalpy is thus H(n) = G(n, t) + T S(n, t).
In closed solutions, the change in enthalpy across a reaction ρ -i.e. the di erence in enthalpy after reaction ρ occurs -is the heat ow from (or to) the thermal reservoir the system is in contact with: (3.29) which can be seen as a statement of the rst law of thermodynamics for closed systems at constant volume [Pippard, 1964].

∆ ρ H(n) = H(n + v ρ ) -H(n) = -Q env ρ ,
Likewise, the change in the entropy of the system across a reaction ρ in a closed system takes the following form:

∆ ρ S(n, t) = S(n + v ρ , t) -S(n, t) = -k B log P (n + v ρ , t) P (n, t) + ∆ ρ H(n) -∆ ρ F (n) T = - ∆ ρ G(n, t) + Q env ρ T , (3.30)
where, in all cases, ∆ ρ stands for the di erence between the state functions after and before the reaction ρ. Finally, the changes in the entropy of the thermal reservoir due to a reaction ρ occurring in a closed system are given by the exchanged heat [Kondepundi andPrigogine, 2015, Rao andEsposito, 2016]:

∆ ρ S env = Q env ρ T . (3.31)
Therefore, taking into account both the entropy of the system and that of the environment, the total entropy production across a reaction is

∆ ρ S tot = ∆ ρ S(n, t) + ∆ ρ S env = - ∆ ρ G(n, t) T . (3.32) 3.3. T CRN 71
A full stochastic trajectory is completely determined by the occurrence of a set of chemical reactions {ρ k } at times {t k }. The trajectory can then be encoded in the instantaneous current for each of the reactions during the dynamics [START_REF] Rao | Conservation laws and work uctuation relations in chemical reaction networks[END_REF]: (3.33) whose probability can be obtained from the CME (3.4).

j ρ (n, τ ) = k δ ρ,ρ k δ n,nt k δ(τ -t k ),
For a trajectory described by Eq. (3.33) -which follows the dynamics of the CME (3.4)entropy production can be written as [Seifert, 2005, Rao andEsposito, 2018]

∆S tot [n t ] = t 0 dτ ∂P (n, τ ) ∂τ + n,ρ (∆ ρ S(n) + ∆ ρ S env )j ρ (n, τ ) = - ∆G[n t ] T , (3.34)
where

∆G[n t ] = t 0 dτ ∂P (n, τ ) ∂τ + n,ρ ∆ ρ G(n, t)j ρ (n, τ ) (3.35)
is the change in free energy along the stochastic trajectory and we use [n t ] to denote a full stochastic trajectory from τ = 0 to time τ = t. The two previous equations are di erences of a state variable along a trajectory that depends on the stochastic occurrences of chemical reactions and the changes in probabilities, these last ones obtained from solving the CME. We emphasise that the previous relations hold only for closed systems, where the entropy production is, indeed, given by the change in free energy, and a generalisation to open non-equilibrium systems will be exposed later.

Finally, it can be proved that the average total entropy di erence along a trajectory is always positive [Seifert, 2005] ∆S tot ≥ 0, (3.36)

which is a statement of the second law of thermodynamics. This also implies that the average stochastic free energy G decreases during a trajectory, as expected in closed systems where no work is being done.

Non-equilibrium thermodynamics of Chemical Reaction Networks

Here we consider a CRN as in the previous section but that has been pushed out of equilibrium by placing it in contact with di erent particle reservoirs, at given chemical potentials. In this case, there are energy uxes owing through the network that require a free-energy expenditure to be maintained, which, as noted in the Preface, is a hallmark of living beings.

Therefore, the reservoirs are exerting a work on the system and a restatement of the rst law of thermodynamics, Eq. (3.29), is needed to account for it. Across a reaction ρ there will be now two di erent terms that contribute to the change in enthalpy [START_REF] Rao | Conservation laws and work uctuation relations in chemical reaction networks[END_REF]]

∆ ρ H(n) = -Q env ρ -Q chem ρ + W chem ρ , (3.37)
where W chem ρ = j µ j v ρ j is the work made by the chemostats across reaction ρ, and Q chem ρ = -T j s j v ρ j is the heat exchanged with the chemostats, with h j -T s j = µ j being the enthalpic C 3. D CRN (h j ) and entropic (s j ) part of the chemostatted chemical potentials µ j , respectively. This is an expression of the rst law of thermodynamics for a system at constant volume on which work is being done by placing it in contact with chemostats at di erent chemical potentials. The decomposition of the heat the system exchanges follows from the identi cation of work performed in the system as W chem ρ = j µ j v ρ j , which is a meaningful and intuitive de nition, particularly if we keep in mind the canonical view of how a biological cell stores energy: in the form of a chemical potential di erence between ATP and ADP.

As before, let us now consider the changes in the entropy of the system across a reaction ρ, which, unlike before, may be driven by the work performed by the chemostats:

∆ ρ S(n, t) = -k B log P (n + v ρ , t) P (n, t) + ∆ ρ H(n) -∆ ρ F (n) T = - ∆ ρ G(n, t) + Q env ρ + Q chem ρ -W chem ρ T . (3.38)
Similarly, the entropy change of the thermal bath is

∆ ρ S env = Q env ρ + Q chem ρ T .
(3.39)

These two relations now include contributions from the chemostats, both in the form of heat and work.

Then, the total entropy production along a stochastic trajectory for a CRN on which chemical work is being exerted is given by [START_REF] Rao | Conservation laws and work uctuation relations in chemical reaction networks[END_REF]]

∆S tot [n t ] = t 0 dτ ∂P (n, τ ) ∂τ + n,ρ (∆ ρ S(n) + ∆ ρ S env )j ρ (n, τ ) = - ∆G[n t ] -W chem [n t ] T (3.40)
where, unless speci ed, di erences in state variables (S tot and G) are taken along a full trajectory and not only across a single reaction and W chem [n t ] is the total chemical work performed along the trajectory. Again, one can prove that the average entropy production is positive [Seifert, 2005, Rao andEsposito, 2018],

∆S tot ≥ 0, (3.41) even in the presence of external work. This yields a bound for the work required to transform a system from one state to another, which reminiscent of the inequalities obtained from equilibrium thermodynamics:

W chem ≥ ∆G . (3.42)
However, ∆G need not be evaluated at thermodynamic equilibrium but can correspond to any non-equilibrium steady state which requires free-energy dissipation in order to be maintained.

A re nement of these ideas can be found in [START_REF] Rao | Conservation laws and work uctuation relations in chemical reaction networks[END_REF], obtained by separating reversible from irreversible work sources. However, for our purposes, the bounds stated in this section are enough and will be applied in the next chapter.

Chapter 4

The role of complex balance in Reaction-Di usion networks

For over a century it has been known that if for each reaction ρ the backward and forward rates are equal, then the steady state of a CRN with MAK satis es the thermodynamic equilibrium condition [Lewis, 1925]. However, the steady state of a CRN with MAK need not be the point of thermodynamic equilibrium, as pointed out rst by R. Wegscheider in 1902 and later, in a more general way, by L. Onsager in 1931 [Onsager, 1931].

Nevertheless, a formalism akin to that of equilibrium thermodynamics can emerge, for MAK but also for more general kinetics, in the sense of being able to de ne a function whose value is minimised by the kinetics of the system and relating its minimum value with the steady state of the system. This is de nitely the case for complex-balanced systems, as it will be shown below.

De nition of complex and complex-balanced network

In the previous chapter, we de ned a complex as the set of chemical species (with their respective particle numbers) that take part in reaction, either as a reactant or as a product. Its most general expression is the vector

z = (z 1 , • • • , z i , • • • ) (4.1)
where the index i runs over all chemical species and z i is an integer coe cient that speci es the number of molecules of the species X i that appear in the complex z. Then, any CRN (as described in Section 3.1) can be represented as a graph whose nodes represent the complexes that take part in the reactions and there is an edge between two complexes if and only if the reaction z m z n exists, where z m and z n denote two di erent complexes (see e.g. In a deterministic CRN, whose kinetics are given by Eq. (3.7), the rate of creation of the complex z is de ned as

J +z = ρ:s ρ =z J +ρ + ρ:r ρ =z J -ρ (4.2)
where the subscript ρ : s ρ = z indicates that the sum is taken over the reactions ρ whose product complex equals the complex z and similarly for ρ : r ρ = z, which denotes that that the sum is 73 C 4. C taken over those reactions whose reactant complex equals the complex z. Analogously, we can de ne the rate of annihilation of the complex z as

J -z = ρ:s ρ =z J -ρ + ρ:r ρ =z J +ρ . (4.3)
Then, a deterministic network is said to have a complex-balanced steady state if its steady state satis es the condition that the creation rate and the annihilation rate of each complex are equal [START_REF] Horn | General mass action kinetics[END_REF]. Mathematically, this means that

J +z = J -z ∀ z.
(4.4)

We have de ned complex-balancing only for a deterministic network but below we will see how this condition will also a ect the stochastic dynamics.

A hierarchy of steady-states

Given the de nitions above, a complex-balanced steady state is only a subset of all the steady states a general CRN can have. However, complex-balanced steady states are still a broad class of steady states which includes, but is not restricted to, detailed-balance steady states.

Therefore, we can order these types of steady states in terms of their generality:

1. The most restrictive condition we can impose to a steady state is that of detailed balance:

J +ρ = J -ρ ∀ ρ, ( 4.5)
condition that corresponds to a system at thermodynamic equilibrium and implies that the rate of the forward reaction equals the rate of the backward reaction, for every reaction ρ.

2. More general than detailed-balanced steady states are complex-balanced steady states, which satisfy that ρ:sρ=z

J +ρ + ρ:rρ=z J -ρ = ρ:sρ=z J -ρ + ρ:rρ=z J +ρ ∀ z, (4.6)
which implies that the creation rate and the annihilation rate of each complex are equal.

In the example of Fig. 4.1, it corresponds to having the net rate (J +ρ -J -ρ ) in reaction 1 4.1. D -75 equal to those of 2 and 3, and, thus, reactions 4 and 5 have to be detailed balanced (zero net rate for each reaction). While it is not the case for the CRN in Fig. 4.1, there are CRN topologies that only accept complex-balanced steady states, as we will see in the Section 4.2.

3. Finally, we have the most general class of steady states, which is just de ned by vanishing time derivatives of the dynamical equation (3.7); condition that, by splitting it into the contributions of each complex z, can be recast into

m z m i   ρ:sρ=zm J +ρ + ρ:rρ=zm J -ρ   = m z m i   ρ:rρ=zm J +ρ + ρ:sρ=zm J -ρ   ∀ i, (4.7)
where m is an index that labels each of the complexes in the network and the integer z m i represents its components. Then, there are no longer constraints between the net rates of each reaction [other than those imposed by Eq. (3.7)]. For Fig. 4.1, this implies that there can appear cycles at steady state where, for example, the species A is created by reaction 1 but annihilated by reaction 4, which breaks complex balance. A generic steady state can allow cycles that cannot be visualised directly from the network representation of the CRN in terms of complexes (such as the graphical representation in Fig. 4.1), while a complex balanced steady state only allows cycles that can be directly visualised from the network representation of a CRN in terms of complexes.

From this hierarchical classi cation, it can be clearly seen that detailed balancing (4.5) implies complex balancing (4.6) and that complex balancing (4.6) implies the system is at steady state (4.7). However, the converse is not true. Indeed, a steady state (4.7) does not necessarily have to be complex balanced (4.6) and a complex-balanced steady state (4.6) is more general than a detailedbalanced one (4.5). Therefore, complex balance is a less restrictive constrain than detailed balance, but it is still less general than a fully generic steady state. Moreover, a complex-balanced network need not be at thermodynamic equilibrium and can sustain net currents at steady state.

Complex balance in networks with mass action kinetics

As a particular case of especial importance, we give the explicit condition for complex-balancing in a network with MAK. The dynamics of a deterministic and ideal CRN [Eqs. (3.7) and (3.8)], with state independent rate constants, can be written as

∂c i ∂t = ρ v ρ i k +ρ i c r ρ i i -k -ρ i c s ρ i i , ( 4.8)
where the non-equilibrium contribution of the chemostats has been absorbed into the rate constants k +ρ and k -ρ . A deterministic and ideal CRN has a complex balanced steady state if for each complex z [START_REF] Anderson | Product-form stationary distributions for de ciency zero chemical reaction networks[END_REF]:

ρ:s ρ =z k +ρ i (c SS i ) r ρ i -k -ρ i (c SS i ) s ρ i = ρ:r ρ =z k +ρ i (c SS i ) r ρ i -k -ρ i (c SS i ) s ρ i , (4.9)
where c SS i is the steady-state concentration of species i 1 . C 4. C

Review of results for complex-balanced systems

In this section we review well-known results from the literature of ideal CRNs, which we will later generalise to non-ideal networks.

It was rst realised half a century ago that CRNs that relax to complex-balanced steady states have a certain number of useful properties. In the seminal work by F. Horn and R. Jackson in 1972, they showed that if a deterministic CRN with MAK accepts a complex-balanced steady state (i.e. Eq. (4.9) is ful lled) then the system is quasi-thermodynamic, which they de ned as a system that minimises a pseudo-Helmholtz function [START_REF] Horn | General mass action kinetics[END_REF]]:

L = i c i log c i -log c SS i -1 , (4.10)
where c SS i is the steady-state concentration of species i. The choice of terminology was due to the fact that Horn and Jackson were looking for CRNs that satis ed a formalism similar to that of equilibrium thermodynamics, where the minimum of a free energy determines the equilibrium point of a system. Indeed, they found the function L, whose minimum speci es the steady state of the system, but this function is not necessarily equal to the equilibrium free energy. This is due to the fact that, as pointed out by Wegscheider and Onsager among others, for a CRN with MAK to relax to equilibrium there are some constraints the rate constants must satisfy, which relate them to the standard-state chemical potentials of the species, as in Eq. (3.15). However, to minimise a pseudo-Helmholtz function like Eq. (4.10), these constraints are not needed and, therefore, it may correspond to a MAK system out of thermodynamic equilibrium.

Shortly thereafter, M. Feinberg showed that, just by examining the topology of the network, one can know if the steady state of the network must be complex balanced [Feinberg, 1972]. Given a CRN, one can de ne the de ciency of the network, δ, which is given by

δ = C -L -S (4.11)
where C is the number of complexes in the network, L is the number of linkage classes or connected components in the network and S is the dimension of the stoichiometric subspace, which amounts to the number of degrees of freedom in the dynamics. If δ = 0, then the system necessarily has a complex-balanced steady state. Therefore, we can assess whether a MAK system will have a pseudo-Helmholtz function as a Lyapunov function only from looking at the topology of the network. Moreover, the de ciency of a network has an interpretation in terms of cycles: it gives the number of cycles that can exist at steady state but cannot be directly visualised from looking at the network graph in terms of the complexes. In other words, it gives the number of emergent or hidden cycles of a CRN [START_REF] Polettini | Dissipation in noisy chemical networks: The role of de ciency[END_REF].

Example 3 To illustrate the concept of de ciency we discuss the CRNs depicted in Fig. 4.2.

First, we compute the de ciency δ of the CRN on the left of Fig. 4.2, which corresponds to a simple enzyme kinetics model with a substrate S, a product P , an enzyme E and an enzyme-substrate complex ES. This CRN has 6 complexes C (including the empty complex ∅ that represents exchanges with particle reservoirs) and two connected components L. Finally, the dimension of the stoichiometric subspace S is 3, The left CRN has de ciency δ = 1 and the right one is zero-de cient δ = 0. Therefore, the network on the right will necessarily have a complex-balanced steady state but this is not the case for the left-hand one.

since the number of enzymes E and enzyme-substrate complexes ES is always conserved. Then de ciency δ = C -L -S = 1, which means that:

1. The steady state of this network need not be complex balanced.

2. There is one emergent cycle that is not represented graphically in the network representation of the complexes (see Fig. 4.2,left CRN). This emergent cycle is: the extraction of a substrate particle S from the particle reservoir, its binding to the enzyme to form a ES complex, conversion to a product P and nal withdrawal of the product P from the system to the particle reservoir.

The second CRN (right CRN of Fig. 4.2) is zero-de cient. It has 5 complexes C, 2 linkage classes L and a stoichiometric subspace, S, of 3 (since there are two conservation laws: D + B + C is constant and the total mass of the system is also constant since there is no exchange with particle reservoirs). Then δ = C -L -S = 0 which means that the steady state of this network is necessarily complex balanced and thus, if the kinetics are MAK, it has a Lyapunov function of the form (4.10). . More recently, it was found by D. F. Anderson, G. Craciun and T. D. Kurtz that, in the stochastic description, the steady-state probability distribution of a complex-balanced CRN with MAK can be written as a product of independent Poisson distributions [START_REF] Anderson | Product-form stationary distributions for de ciency zero chemical reaction networks[END_REF]. Under this constraint, the steady-state distribution π CB of the network with dynamics (3.4) and propensity functions (3.5) is given by

π CB (n) = M i=1 (c SS i ) n i n i ! e -c SS i , (4.12)
where c SS i are the steady-state solutions of the complex-balanced networks in the deterministic limit [this is, Eqs. (3.7) and (3.8)]. Note that the expression (4.12) for the steady-state probability distribution is valid only for vectors n satisfying all the conservation laws the CRN may have, otherwise, its value is zero. This distribution includes the equilibrium distribution of the CRN (the Boltzmann distribution for an ideal system), but is not limited by thermodynamic equilibrium and may be applied to cases far from equilibrium. Furthermore, they also proved a more general C 4. C result for propensity functions of the form

f +ρ (n) = k +ρ θ(n) θ(n -r ρ ) , (4.13)
which, in complex-balanced systems, yields a steady-state distribution with the following form [START_REF] Anderson | Product-form stationary distributions for de ciency zero chemical reaction networks[END_REF]:

π CB (n) = M θ(n) i=1 (c SS i ) n i , (4.14)
where M is a normalisation constant and θ is a function that maps the vector of integer numbers n into a single real-valued positive number.

To conclude with this review of the literature, recently D. F. Anderson and T. D. Nguyen genralised the Lyapunov function of complex-balanced MAK systems [Eq. (4.10)] to networks with product-form stationary distributions [START_REF] Anderson | Results on stochastic reaction networks with non-mass action kinetics[END_REF]. In the following, we will further generalise some of these results and will give an interpretation of them in terms of thermodynamic quantities.

Complex-balanced steady-state distributions

Closely inspecting Eqs. (4.13) and (4.14) and identifying θ(n) with the Boltzmann factor e -βF (n) , one can already get a hint of the connection with the thermodynamical constraints imposed on the propensity functions f ρ in the previous chapter. In this section, we will look at a slight generalisation of the probability distribution (4.14), by considering a complex-balanced CRN, with propensity functions given by Eq. (3.25). However, in the following, we will assume that the function g ρ in Eq. (3.25) is the same for all reactions, that is, g ρ = g. Therefore, we do not require that the function g is a constant, but we do need the fact that the rate constants of all reactions depend on the state n in the same way, which implies that g for all reaction must be the same function. This simpli cation allows us to prove the result and the consequences of relaxing this hypothesis will be examined in the following Chapter 5.

Under these conditions, we claim that, for complex-balanced CRNs the steady-state solution of the stochastic dynamics (3.4) with the propensity functions (3.25) takes the form:

π CB (n) = e -β[F (n)+ i μi n i ] Z , ( 4.15)
where μi are parameters that will depend on the chemostats to which the system is connected (and that push the system out of equilibrium) and the reaction constants of the CRN. In fact, the parameter μi is related to c SS i in the previous section, but we have changed its form to better illustrate its thermodynamic meaning. Importantly, these μi parameters can be obtained from the ideal and deterministic CRN, greatly simplifying the task of obtaining analytically the steadystate distribution of the system. By analogy with equilibrium statistical physics, we denote by Z the normalisation factor of the probability distribution.

The strategy we will follow to prove the result is similar to that used by Anderson, Craciun and Kurtz in Ref. [START_REF] Anderson | Product-form stationary distributions for de ciency zero chemical reaction networks[END_REF], where the key point is to simplify the complex-balanced 4.3. C --79 steady-state relation until we obtain that of a deterministic CRN modelled with MAK, given in Eq. (4.9). Then, we will be in a position to compare our parameters μi to the steady-state concentrations of the chemical species in the deterministic MAK case and write down the full probability distribution at steady state for complex-balanced non-ideal networks.

Proof

The CME (3.4), with rates of the form (3.25) and g ρ = g for all reactions ρ, at steady state becomes

ρ kρ g(n -s ρ )e β[F (n-vρ)-F (n-sρ)+ j r ρ j µ j ] P (n -v ρ ) + ρ kρ g(n -r ρ )e β[F (n+vρ)-F (n-rρ)+ j s ρ j µ j ] P (n + v ρ ) = ρ kρ g(n -r ρ )e β[F (n)-F (n-rρ)+ j r ρ j µ j ] + kρ g(n -s ρ )e β[F (n)-F (n-sρ)+ j s ρ j µ j ] P (n).
(4.16)

Dividing the previous expression by P (n) and inserting an ansatz for the solution of the form (4.15) we obtain:

ρ kρ g(n -s ρ )e β[F (n)-F (n-sρ)+ i v ρ i μi + j r ρ j µ j ] + g(n -r ρ )e β[F (n)-F (n-rρ)-i v ρ i μi + j s ρ j µ j ] = ρ kρ g(n -r ρ )e β[F (n)-F (n-rρ)+ j r ρ j µ j ] + g(n -s ρ )e β[F (n)-F (n-sρ)+ j s ρ j µ j ] . (4.17)
We can write this relation in terms of a summation over each of the complexes z separately

z ρ:sρ=z kρ g(n -s ρ )e β[F (n)-F (n-sρ)+ i v ρ i μi + j r ρ j µ j ] + z ρ:rρ=z kρ g(n -r ρ )e β[F (n)-F (n-rρ)-i v ρ i μi + j s ρ j µ j ] = z    ρ:rρ=z kρ g(n -r ρ )e β[F (n)-F (n-rρ)+ j r ρ j µ j ] + ρ:sρ=z kρ g(n -s ρ )e β[F (n)-F (n-sρ)+ j s ρ j µ j ]    , (4.18)
where the subscript ρ : s ρ = z denotes that the sum runs only over all reactions ρ whose product complex s ρ equals the complex z. Until now we are considering an arbitrary steady state.

Introducing the complex balance constraint, we can impose that, separately, the uxes across all complexes vanish (that is J +z = J -z ∀ z), which, can be rewritten in the following form, given that the complex z is xed:

g(n -z)e β[F (n)-F (n-z)]    ρ:sρ=z kρ e β[ i v ρ i μi + j r ρ j µ j ] + ρ:rρ=z kρ e β[-i v ρ i μi + j s ρ j µ j ]    = g(n -z)e β[F (n)-F (n-z)]    ρ:rρ=z kρ e β j r ρ j µ j + ρ:sρ=z kρ e β j s ρ j µ j    . (4.19) C 4. C
Then, Eq. (4.18) will be satis ed if for each complex, z, Eq. (4.19) is sati ed. We can now divide both sides by g

(n -z) exp(β[F (n) -F (n -z)]) to obtain ρ:sρ=z kρ e β[ i (z i -r ρ i )μ i + j r ρ j µ j ] + ρ:rρ=z kρ e β[-i (s ρ i -z i )μ i + j s ρ j µ j ] = ρ:rρ=z kρ e β j r ρ j µ j + ρ:sρ=z kρ e β j s ρ j µ j , (4.20)
where we have substituted v ρ = s ρ -r ρ and, according to each of the sums, one of this terms can be replaced by the complex z. Finally, given that z is xed in the previous equation, we can divide the whole expression by exp(β i z i μi ), yielding

ρ:sρ=z kρ e β[-i r ρ i μi + j r ρ j µ j ] + ρ:rρ=z kρ e β[-i s ρ i μi + j s ρ j µ j ] = ρ:rρ=z kρ e β[ j r ρ j µ j -i z ρ i μi ] + ρ:sρ=z kρ e β[ j s ρ j µ j -i z ρ i μi ] . (4.21)
At this point, if we make the following de nitions:

c SS i = exp(-β(μ i + µ 0 i )), (4.22
)

k +ρ = kρ exp(β[ i r ρ i µ 0 i + j r ρ j µ j ]), (4.23)
k -ρ = kρ exp(β[ i s ρ i µ 0 i + j s ρ j µ j ]), (4.24)
then we recover the complex balance condition (4.9) for an ideal and deterministic network. Hence, a CRN for which the deterministic steady-state is complex balanced, accepts a steadystate solution of the form (4.15) for its stochastic and non-ideal version. The parameters μi that appear in the steady-state distribution can be obtained using Eq. (4.22) and solving the ideal and deterministic CRN for the steady state concentrations (c SS i ), with rates constants as speci ed in Eqs. (4.23) and (4.24), where the standard-state chemical potentials and non-equilibrium contributions from the chemostats have been absorbed into the rate constants, as is customary in MAK.

Therefore, the steady-state distribution of a non-ideal complex-balanced CRN has the form of an e ective Boltzmann distribution, where the standard-state chemical potentials µ 0 i are shifted by an amount μi . Physically, it is interesting the fact that, for complex-balanced steady states, the free-energetic contribution decouples from the non-equilibrium terms, embodied in the shifted chemical potentials μ. This may not be true for more general steady-states. This result is similar to Theorem 6.6 of Ref. [START_REF] Anderson | Product-form stationary distributions for de ciency zero chemical reaction networks[END_REF], whose main result is the probability distribution (4.14), but here we have generalised it slightly to include rates of the form (3.25), which includes the function g that could be of interest in phase-separated systems as it modulates the rates depending on the environment. Moreover, our approach clari es the thermodynamic origin of the rates (3.25) and the steady-state distribution (4.15), as is shown in the next subsection. Example 4 Let us consider the following CRN (see Fig. 4.3 for a graphical representation):

A + D B B C C A + D, (4.25)
with a free energy taken from a regular solution theory (see Appendix C.1). Given that the number of particles is not conserved, for simplicity, we assume that the solvent volume is conserved and allow the total volume to vary (alternatively, we could have chosen A and D particles as occupying half the volume as C or B). Thus, V = N + n sol + n D , where N is the number of conserved particles N = n A + n B + n C and n sol is the number of solvent particles. Moreover, given any initial condition, the variable n D will always be ensalved to the variable n A as they are both created and annihilated together, keeping constant the initial particle di erence between the two. For simplicity, we assume that the system is driven out of equilibrium solely by imposing a non-equilibrium chemical potential di erence ∆µ neq = log(5/2)β -1 in the transition from C to A + D, which is equivalent to assuming that the reaction C A + D actually is C + E A + D + F , where E and F are chemostatted species with chemical potentials satisfying

µ E -µ F = log(5/2)β -1 .
We assume all the reaction constants are equal to each other and note that the network is necessarily complex-balanced: unimolecular reactions are always complex balanced and here particles A and D are only created or destroyed via the complex A + D, therefore the steady-state condition yields that the creation and annihilation rates of the complex A + D must be equal. Equivalently, one could compute the de ciency of the network and check that δ = 0, hence the steady state must be complex balanced. Then, the steady-state of the CRN (4.25), modelled stochastically with dynamics (3.4) and propensity functions (3.25), can be obtained from the following deterministic rate equations (that include entropic terms and the non-equilibrium driving, but not any other term of the free energy):

dn A dt = n B + 2.5n C -2n A n D dn B dt = n C + n A n D -2n B dn C dt = n B + n A n D -3.5n C .
π CB (n) = e -β[F (n)+ i μi n i ] Z , ( 4.26)
once the free energy F has been speci ed and the normalisation constant Z computed. Note that there are two conservation laws (N = 40 and n A -n D = 5), and four chemical species, hence, in practice, π neq (n) is a distribution with only two independent variables.

Let us consider the following regular solution free energy, that comprises standard-state chemical potentials (taken with respect to that of species B), interactions (parametrised by χ) and entropic terms:

F (n A , n B , n C , n D ) = β -1 [log(n A !n B !n C !n D !n sol !) -log V !] + n A µ 0 A + n C µ 0 C + n D µ 0 D + χ n A n C V , (4.27)
where V = n A + n B + n C + n D + n sol and n sol is the number of solvent molecules. For the parameters χ = 10β -1 and µ

0 A = µ 0 C = µ 0 D = log(2)β -1
, we obtain the steady-state probability depicted in Fig. 4.4, that matches closely the one obtained from a simulation of the same CRN using the Gillespie algorithm [Gillespie, 1977]. Simulations were started in parallel from random poissonian initial conditions (provided they satisfy the constraints above) and the samples were obtained after allowing the simulations to relax to the steady state.

Thermodynamic interpretation of μ

The direct consequence of the non-equilibrium driving of the complex-balanced CRN is the appearance of a term μ in the stationary distribution (4.15). E ectively, this term shifts the standardstate chemical potential µ 0 by a constant. However, the thermodynamic implications of this shift in the standard-state chemical potentials can be further clari ed.

L -83

As shown in Section 3.3.2, one can bound the work required to make a non-equilibrium transformation of the system by the di erence in free energies between these non-equilibrium states. Here, we consider a simpler case of obtaining a non-equilibrium complex-balanced steady state from an initial equilibrium state, by connecting the system to several chemostats. Then, in the stochastic description, the average work required to do that is bounded by

W chem ≥ G f -G eq , (4.28)
where G is the stochastic free energy de ned by Eq. (3.27).

At equilibrium, the probability distribution is just the canonical distribution (or grand-canonical distribution in case the particle numbers are not conserved), which is, up to a normalisation constant (the partition function Z), the exponential of the free energy F . Thus, from Eq. (3.27), G eq = -k B T log Z, as expected from equilibrium statistical physics and where Z is the equilibrium partition function. However, for a non-equilibrium complex-balanced state,

G CB = n π CB (n) i μi n i -k B T log Z CB = i μi n i -k B T log Z CB , (4.29)
where Z CB is the complex-balanced normalisation constant, which is not necessarily equal to the equilibrium partition function Z, and n i is the average number of i particles, according to the complex-balanced distribution π CB .

Therefore, the work required for such a transformation is bounded by

W chem ≥ i μi n i f + k B T log Z Z CB , (4.30)
where n i f is the average number of particles i at the nal non-equilibrium complex-balanced steady state. In the previous relation, the ratio of partition functions is reminiscent of the equilibrium bound for a process, which is only attained if the transformation is fully reversible (∆S tot = 0) during the transformation. However, the appearance of the i μi n i f is exclusively due to the non-equilibrium work performed on the system, since μi = 0 in the absence of such work.

Thus, this inequality suggests a thermodynamic interpretation of the shifted chemical potentials μ, since the bound depends crucially on them and is the only inherently non-equilibrium contribution to the bound.

Lyapunov function for complex-balanced steady states

Under fairly general conditions, the logarithm of the steady-state probability distribution in the stochastic model is a Lyapunov function of the deterministic one, provided an appropriate scaling is used to bridge between the di erent volume and particle-number scales [START_REF] Ge | Mesoscopic kinetic basis of macroscopic chemical thermodynamics: A mathematical theory[END_REF]. In particular, this has been shown for ideal and complex-balanced CRNs [START_REF] Anderson | Lyapunov functions, stationary distributions, and non-equilibrium potential for reaction networks[END_REF]. Here, we show that in non-ideal complex-balanced CRNs the following function is minimised by the dynamics

L(c) = -lim V →∞ 1 V log(π CB (n)) = β f (c) + i μi c i + log(Z)/V, (4.31) C 4. C
where c i = n i /V and f (c) = F (n)/V . Indeed, this approach generalises to a certain extent the results of [START_REF] Anderson | Results on stochastic reaction networks with non-mass action kinetics[END_REF] for product-form stationary states of CRNs. We abuse slightly the terminology by calling this function a Lyapunov function as we do not prove it is strictly positive or bounded from below, we will only prove that it is minimised by the dynamics. Nevertheless, the fact that this function is minimised is su cient for our purposes.

Since the normalisation constant Z does not depend on time (it only depends on the nonequilibrium steady state), the time derivative of the Lyapunov function can be written as

dL dt = i ∂L ∂c i ∂c i ∂t =β i (µ i + μi ) ρ v ρ i k ρ g(c)(e β[ i r ρ i µ i + j r ρ j µ j ] -e β[ i s ρ i µ i + j s ρ j µ j ] ) , (4.32)
where we have used Eq. (3.7) with currents given by Eq. (3.26). Note we have used two di erent dummy indices, i and i , to make it clear to which sum each variable corresponds to. We can recast the previous expression into

dL dt = ρ i k ρ g(c)(µ i + μi )(s ρ i -r ρ i )e β[ i r ρ i (µ i +μ i )-i r ρ i μi + j r ρ j µ j ] + ρ i k ρ g(c)(µ i + μi )(r ρ i -s ρ i )e β[ i s ρ i (µ i +μ i )-i s ρ i μi + j s ρ j µ j ] , (4.33)
by adding and subtracting the μ-terms in the arguments of the exponential functions. We now can repeatedly apply the inequality 2 e a (b -a) ≤ e b -e a to the sums of chemical potentials to arrive at:

dL dt ≤ ρ k ρ g(c)e β[ j r ρ j µ j -i r ρ i μi ] e β i (µ i +μ i )s ρ i -e β i r ρ i (µ i +μ i ) + ρ k ρ g(c)e β[ j s ρ j µ j -i s ρ i μi ] e β i (µ i +μ i )r ρ i -e β i s ρ i (µ i +μ i ) . (4.34)
This previous expression can now be separated in terms of the di erent complexes in the system:

dL dt ≤ z g(c) ρ:s ρ =z k ρ e β[ j r ρ j µ j -i r ρ i μi + i (µ i +μ i )s ρ i ] - ρ:r ρ =z k ρ e β[ j r ρ j µ j -i r ρ i μi + i r ρ i (µ i +μ i )] + ρ:r ρ =z k ρ e β[ j s ρ j µ j -i s ρ i μi + i (µ i +μ i )r ρ i ] - ρ:s ρ =z k ρ e β[ j s ρ j µ j -i s ρ i μi + i s ρ i (µ i +μ i )] ,
(4.35)

where, invoking complex balancing as in the previous section and dividing by exp[β i z i (µ i +

2 The inequality stems from 1 + x ≤ e x , ∀ x ∈ R and setting x = b -a. We use it by identifying (for the rst term) 4.5. L 85 μi )], we have an equality for each complex z

a = i (µi + μi)r
ρ:s ρ =z k ρ e β[ j r ρ j µ j -i r ρ i μi ] - ρ:r ρ =z k ρ e β[ j r ρ j µ j -i r ρ i μi ] + ρ:r ρ =z k ρ e β[ j s ρ j µ j -i s ρ i μi ] - ρ:s ρ =z k ρ e β[ j s ρ j µ j -i s ρ i μi ] = 0, (4.36)
equality that is given by the complex-balance condition for deterministic MAK networks [Eq. (4.9)] together with Eqs. (4.22) to (4.24). Summing over all complexes yields the desired inequality: dL dt ≤ 0, (4.37) and, thus, L is minimised (or remains unchanged) along a trajectory for a deterministic complexbalanced non-ideal CRN.

Therefore, unlike in Section 2.2.5 where we minimised the free energy F since we were looking for the equilibrium state, here we minimise the function L that includes the free energy density f and the standard-state chemical potentials µ 0 i have been shifted by μi due to the work performed on the system by the chemostats. Consequently, the CRN does not relax to thermodynamic equilibrium and, at steady state, is able to sustain chemical currents.

Lyapunov functionals for spatially heterogeneous solutions

To be able to describe phase-separating systems we need to take into account the spatial degrees of freedom. In the deterministic description, concentrations are now a function of space, c i (x), within a volume Ω and the free energy is a functional of these concentrations, F [c(x)].

The time derivative of the concentration elds c i (x) is given by the following reaction-di usion equation

dc i (x) dt = -∇ • J i + ρ v ρ i (J +ρ -J -ρ ), (4.38)
where the rst term in the right-hand side (RHS) of the equation corresponds to di usion and the second one to the chemical reactions. As in the linear irreversible thermodynamics framework [Groot andMazur, 1983, Kondepundi andPrigogine, 2015], the driving force of the di usion current J i is the gradient of chemical potentials ∇µ i . More precisely, the di usion currents take the following form J i =j M ij ∇µ j , where M ij is the mobility matrix. We assume noux boundary conditions J i | x=∂Ω = 0 for the non-chemostatted species, where ∂Ω denotes the boundaries of the volume Ω.

We now proceed to show that the Lyapunov function obtained above for a complex-balanced CRN is also a Lyapunov functional if we consider the di usion of solutes and solvent within the solution. Then, the time evolution of the functional L(x)dx becomes

dL dt dx = i dx dL dc i dc i dt = i dx β(µ i (x) + μi ) -∇ • J i + ρ v ρ i (J +ρ -J -ρ ) , (4.39) C 4. C
where µ i (x) = δF [c(x)]/δc i (x) is the local chemical potential and, in the second equality, Eq. (4.38) has been used. From the results in the previous section, we know that the contribution of the second term in the square brackets of the RHS of the equation to the time derivative of L is negative, since it corresponds to the chemical reactions. Therefore, to prove that L(x)dx is minimised by the dynamics is su cient to show that the contribution of the rst term in the RHS is also negative.

We note that the contribution of di usion to the time evolution of L can be written as

i dx (µ i (x) + μi )∇ • J i = i dx [∇ • {(µ i (x) + μi )J i } -∇(µ i (x) + μi ) • J i ] . (4.40)
For the rst term, we have that

dx∇ • {(µ i (x) + μi )J i } = ∂Ω dS(µ i (x) + μi )J i • n = 0, (4.41)
where n is the verctor normal to the surface ∂Ω. In the rst equality we have used the divergence theorem and, in the second one, the fact that we are considering no-ux boundary conditions, although it may be generalised to other appropriate boundary conditions. Regarding the second term in Eq. (4.40), we have that

-∇(µ i (x) + μi ) • J i = -∇µ i (x) • J i = Ṡdi ≥ 0, (4.42)
where Ṡdi is the entropy production rate due to di usion. In the rst equality we have used the fact that μi are constants (hence ∇μ i = 0). Then, we can identify the remaining terms as the entropy production due to the di usion process, which necessarily has to be greater or equal to zero, provided the Onsager reciprocal relations for the mobility matrix M ij are met [Groot andMazur, 1983, Kondepundi andPrigogine, 2015].

We nally obtain that dL dt dx ≤ 0, (4.43)

as both the contribution from chemical reactions and from di usion are negative. Thus, a nonideal complex-balanced reaction-di usion system will minimise the functional Ldx, which we can now use to obtain information about the steady state, as done with the free energy F for systems that relax to thermodynamic equilibrium.

Example 5 Let us now consider the following CRN (see Fig. 4.3 for a graphical representation):

A B B C C A, (4.44)
with a free energy taken from a regular solution theory, as before. Again, for simplicity, we assume that the system is driven out of equilibrium solely by imposing a non-equilibrium chemical potential di erence ∆µ for the transition from C to A.

We take all the reaction constants equal to each other and note that the network is necessarily complexbalanced because all chemical reactions are unimolecular. Then, the steady-state of the CRN (4.44), modelled stochastically with dynamics (3.4) and propensity functions (3.25), can be obtained from the following 4.6. P -87 deterministic rate equations (that include entropic terms and the non-equilibrium driving, but not any other term of the free energy):

dc A dt = c B + e -β∆µ c C -2c A dc B dt = c C + c A -2c B dc C dt = c B + c A -(1 + e -β∆µ )c C .
The solution of the above system at steady state is

c SS A = c SS C 1 2 + e -β∆µ 2 3 c SS B = c SS C 1 2 + e -β∆µ 1 3 + 1 2 .
Identifying, c SS i = e -β μi , and noting that we can express them with respect to that of species C, we get the values of shifted non-equilibrium chemical potentials μi and we can obtain the Lyapunov function of the system

L(c) = βf (c) -c A log 1 + 2e -β∆µ 1 3 -c B log 1 2 + 1 2 + e -β∆µ 1 3 -log Z, (4.45)
where the non-equilibrium partition function has not been explicitly computed but, since it is a constant term along the dynamics, it will not alter the location of the minima of L in the space of concentrations c.

Let us consider the following regular solution free energy in the deterministic limit, that comprises interaction and entropic terms and includes Cahn-Hilliard terms in order to model smoothly the potentially non-uniform system:

f (c A , c B , c C ) = β -1 α=A,B,C c α log c α + χ 22 c 2 A + χ 12 c A c B + i,i κ i,i (∇c i ∇c i ), (4.46)
where χ 22 = -2, χ 12 = -7 and κ i,i are the Cahn Hilliard coe cients [START_REF] Cahn | Free energy of a nonuniform system. I. Interfacial free energy[END_REF], which, as introduced in Section 2.2.1, set the free-energetic cost of interfaces and their width. This free energy and the reaction-di usion equation (4.38) can describe the dynamics of the system, which minimise its Lyapunov functional (4.45), as seen in the Figure 4.5 for a one dimensional system. In particular, the upper panel shows that the Lyapunov functional (4.45) is minimised by the numerical time evolution of the reaction-di usion system, the middle panel shows the concentration pro les at steady state and the lower panel the net chemical current from species C to A at steady state as a function of the position x.

Phenomenology of a complex-balanced mixture

Since the Lyapunov functional for complex-balanced systems found in the previous section is minimised by the dynamics, it carries plenty of information about the steady-state in which the macroscopic system will settle. The phase diagram of a mixture, in the case of phase separation, states whether for a certain parameter set the steady state of the system is homogeneous or there are multiple di erent phases. Thus, the phase diagram follows from the Lyapunov functional C 4. C (4.38) forward in time until it reaches a steady state. Top panel: evolution of the value of the Lyapunov functional with time. This functional is minimised with the dynamics but it does not attain the value 0 for two reasons: First, the constant termlog(Z) has been neglected. Second, the dynamics may not have reached the global minimum of L but only a local one. Middle panel: Steady-state con guration of the system, displaying the volume fractions of each species as a function of the spatial coordinate φ(x). In this case, the system reaches at a steady state where there is phase separation. Lower panel: Net reaction ux at steady state for the third reaction in the CRN (4.44) as a function of the spatial coordinate, showing that the system is not detailed-balanced. Note that, while in the upper panel, the horizontal dimension refers to time, it is no longer the case for the lower ones, where it represents the space in a one-dimensional system.

4.6. P -89 (4.31) in the same way it would do for a classical free energy: the phase coexistence lines are obtained from the minimisation of L with the appropriate constraints (particle conservation, in most of the cases we consider).

Phase diagram of a non-ideal complex-balanced solution

In this section we will consider a non-ideal solution with the CRN (4.44) in the deterministic limit, and obtain its phase diagram. Therefore, we need to minimise the Lyapunov functional found in Example 5 subject to the particle conservation constraint

V φ N = dx[φ A (x) + φ B (x) + φ C (x)],
where φ i (x) are the continuous volume fraction elds of each of the species at point x in space (as is customary when discussing phase separation) and φ N is the overall volume fraction of solutes.

Then, the Lagrangian that needs to be minimised is

L = dx L(φ) -λ V φ N -dx [φ A (x) + φ B (x) + φ C (x)] , (4.47)
where λ is the Lagrange multiplier associated with the conservation of solute particles. If we assume the bulk free energy of the homogeneous phases is large compared to the free-energetic cost of the interfaces, then we can drop the continuity of the concentration elds in space and have a simpler description of the system only in terms of two homogeneous but di erent phases. Considering a nite interface would, in general, yield a small correction to the phase diagram that we nd below. Then, we need to minimise the following Lagrangian:

L = V (1) L(φ (1) ) + (V -V (1) )L(φ (2) ) -λ V φ N -V (1) (φ (1)
A + φ

(1)

B + φ (1) C ) -(V -V (1) )(φ (2)
A + φ

(2)

B + φ (2) C ) , (4.48)
where the superscripts (1) and (2) stand for the di erent phases, V (1) is the volume of one of the phases, and we recall that L is the one obtained in Example 5. For simplicity, we consider the following free energy

f (φ) = α φ α log φ α -χφ 2 A . (4.49)
Minimisation of the Lagrangian (4.48), yields the phase diagram in Fig. 4.6 (for details on the minimisation, see Appendix C.3). Regions 2 and 3 in the phase diagrams correspond to parts of the parameter space where phase separation occurs, as it can be seen from the concentration pro les of the upper panels, while in regions 1 and 4 the homogeneous con guration is stable. These results are reminiscent of equilibrium phase diagrams but we recall that here there are net chemical uxes being maintained at steady state due to the non-equilibrium driving. In fact, from the form of the free energy (4.49), we see that species A drives phase separation. Thus, since the overall concentrations are modulated by the non-equilibrium work ∆µ, whether the steady state displays one phase or a coexistence of phases also depends on the value of ∆µ.

General features of the steady state

As seen above, the phenomenology of a complex-balanced system does not change much with respect to a non-ideal solution in equilibrium, but the non-equilibrium terms can move the system Upper panels: Each of the panels has associated a number which states their corresponding point in the phase diagrams. These panels depict the steadystate volume fraction pro les of the system in a one dimensional space, x, from 0 to L. These results were obtained by numerically integrating in time the reaction-di usion equations, as in Fig. 4.5, with the parameters set by their point in the phase diagram. Note that, in the upper panels, the width of the interface is nite to ease the numerical integration, unlike in the lower panels, where the phase diagram was obtained for vanishing interfaces.

4.6. P -91 across the phase diagram, enabling the switching on and o of phase separation. This can be clearly seen from the minimisation of L [Eq. (4.47)] that involves taking functional derivatives with respect to the volume fractions φ i (x): (4.50) which implies that the chemical potential is homogeneous (since μi and λ are constants) as it would be the case for a mixture at thermal equilibrium [cf. Section 2.2.5, and, in particular Eq.

δL δφ i (x) = µ i (x) + μi -λ = 0,
(2.38)]. Hence, if we de ne non-equilibrium patterning as the type of steady-state patterning characterised by displaying di usion currents, it is excluded for solutions driven out of equilibrium by a complex-balanced CRN since there cannot be any di usion currents at steady state (the driving force of the di usion currents ∇µ i = 0). Nevertheless, as shown in Example 5, the system can sustain local chemical currents at steady state, due to the chemical work supplied by the chemostats. However, these chemical currents do not couple to the di usion process to produce di usion currents at steady state.

One of the consequences of our analysis is the fact that breaking complex-balance is a minimal requirement for non-equilibrium pattern formation in non-ideal solutions with local interactions, at least when modelled in a thermodynamically consistent way, which generalises previous considerations [Carati andLefever, 1997, Kirschbaum andZwicker, 2021]. Nevertheless, patterning in non-ideal mixtures can be achieved by the addition of surfactants or long-range interactions, which can yield states of suppressed Ostwald ripening [Tarzia andConiglio, 2006, Pham et al., 2010].

Our analysis highlights the fact that in non-ideal media complex balance (as de ned in the previous chapters) can be broken in two di erent ways. Breaking complex balance can be achieved by a suitable CRN topology, e.g. [START_REF] Carati | Chemical freezing of phase separation in immiscible binary mixtures[END_REF], but, in a system that exhibits coexistence between two di erent phases, one can also break complex balance by allowing the rates of the reactions to depend on their local environment in di erent ways (for example, in Ref. [START_REF] Kirschbaum | Controlling biomolecular condensates via chemical reactions[END_REF] they found a patterned steady-state by allowing one reaction to depend on the concentration of an enzyme that localises preferentially at one of the phases). Here, to prove the form of the steady-state distribution and the Lyapunov function for complexbalanced systems, we assumed that, even if the reaction were allowed to depend on their environment via a function g ρ , this function had to be the same for all reactions. Hence, if di erent reactions have di erent g ρ -functions, our results may not hold and, as demonstrated in Ref. [START_REF] Kirschbaum | Controlling biomolecular condensates via chemical reactions[END_REF], patterning may occur.

In biological cells, phase separation has been hypothesised to perform many functions, such as, accelerating biochemical reactions within the condensate irrespective of the rate of the reaction in the dilute phase [START_REF] Lyon | A framework for understanding the functions of biomolecular condensates across scales[END_REF]. Our work also implies that, in order to accelerate a reaction only in one of the phases, breaking complex balance is necessary. Indeed, in a complex-balanced system, the chemical potential of every species is the same in both phases (in the condensate and outside it and, in fact, according to Eq. (4.50), it is perfectly homogeneous across space). Then, given that the force driving the chemical reactions are the chemical potentials, the reaction rates in both phases are related, making it impossible to regulate the rates of chemical reactions in each phase independently and suggesting that breaking complex balance in one of the two ways outlined C 4. C above is crucial. However, there are other functions that phase separation in biological cells can perform with complex-balanced CRNs or even at thermodynamic equilibrium, like sequestration of molecules or noise bu ering [START_REF] Klosin | Phase separation provides a mechanism to reduce noise in cells[END_REF].

Complex balance is known to be a key feature of the CRN that determines not only its behaviour [START_REF] Anderson | Product-form stationary distributions for de ciency zero chemical reaction networks[END_REF] but also its thermodynamic properties [START_REF] Polettini | Dissipation in noisy chemical networks: The role of de ciency[END_REF]. Here, we have further stressed the connection between the topology of the CRN and the thermodynamically consistent structure of the physical system to generalise results from ideal CRNs and explore non-equilibrium thermodynamics of complex-balanced non-ideal CRNs. However, for non-complex-balanced systems little is known still and, given our results, further research regarding the behaviour of this type of networks will be of the utmost importance both from the physical and the biological point of view. Some steps in this direction will be taken in the chapter that follows.

Chapter 5

Beyond complex balancing?

In this work we have shown that for chemically-reactive non-ideal solutions we can obtain results analogous to those of ideal CRNs, provided the system is modelled in a thermodynamically consistent way. This implies that the rates of the chemical reactions will be a ected by the interactions between the species in the system and, therefore, MAK is no longer a valid description for the dynamics of the CRN. In the previous chapter, generalising MAK for a non-ideal solution, we obtained the steady-state probability distribution for a stochastic complex-balanced CRN and the Lyapunov function of its deterministic counterpart, which speci es the phase diagram of the system.

Our results are of particular importance for non-equilibrium phase-separating systems, where a wide variety of non-equilibrium behaviour can arise [START_REF] Weber | Physics of active emulsions[END_REF]. However, they were limited to a class of CRNs known as complex-balanced and, as discussed above, many of the inherently non-equilibrium phenomena can only occur for systems that are not complex balanced. Examples of these are: the modulations of reaction rates in one phase (e.g. a phase-separated condensate in a cell) independently of the rest of the system [START_REF] Lyon | A framework for understanding the functions of biomolecular condensates across scales[END_REF] or non-equilibrium patterning [START_REF] Carati | Chemical freezing of phase separation in immiscible binary mixtures[END_REF]]. Thus, going beyond complex-balanced networks seems imperative to understand the phenomenology non-ideal solutions may present in biological cells. This is a very challenging task and, in the present chapter, we will limit ourselves to analysing a particular example and provide directions towards which future research should be directed.

Patterning when complex balance is broken

As an example of behaviour that can arise when complex balance is broken, we consider a simplied version of the model analysed in Ref. [START_REF] Kirschbaum | Controlling biomolecular condensates via chemical reactions[END_REF]: A ternary mixture (two chemical species and a solvent) where the two chemical species are exchanged via chemical reactions of the form A B.

First we will consider the equilibrium behaviour of the system. Afterwards, we will compare it to its behaviour out of thermodynamic equilibrium while restricted to complex-balanced networks. Finally, we will explore the e ects of breaking complex balance out of equilibrium.

Equilibrium behaviour of the mixture

Let us assume that the system is not chemically reactive and is well described by the following equilibrium free energy density

f k B T = c A log(c A ) -c A + c B log(c B ) -c B -χc 2 A + κ 2 |∇c A | 2 + |∇c B | 2 , (5.1)
where κ are the Cahn-Hilliard coe cients that penalise the creation of interfaces. The eigenvalues of the Hessian matrix of the free energy density yield the stability of the homogeneous solution and can be used to obtain the parameter region for spinodal decomposition (see Chapter 2, Section 2.3.1). Considering a homogeneous system (where ∇c i = 0), we obtain that the smallest (most negative) eigenvalue of the Hessian matrix is given by

λ -= 1 c A -2χ, (5.2)
which controls the stability of the system and, hence, whether spinodal decomposition occurs. Therefore, when 2χ > 1 c A , the homogeneous con guration of the system is unstable and spinodal decomposition will occur.

Complex-balanced behaviour

Let us now consider the case of the exchange A B that can happen via two di erent reaction pathways: one with a non-equilibrium contribution from the chemostats ∆µ and the other one without any external work. Then, one can write the resulting reaction di usion-equations for the system:

∂c A ∂t =D A ∇ 2 µ A + k 1 (e β(µ B +∆µ) -e βµ A ) + k 2 (e βµ B -e βµ A ) (5.3) ∂c B ∂t =D B ∇ 2 µ B -k 1 (e β(µ B +∆µ) -e βµ A ) -k 2 (e βµ B -e βµ A ), (5.4)
where µ A = δf /δc A and µ B = δf /δc B are the chemical potentials of species A and B, respectively. A concentration homogeneous pro le is a steady-state solution to Eqs. (5.3) and (5.4), provided the chemical potentials satisfy

e βµ A = k 1 e β∆µ + k 2 k 1 + k 2 e βµ B .
(5.5)

In order to analyse the stability of this solution we perform a Linear Stability Analysis (LSA) (for an introduction to LSA see e.g. [Kondepundi and Prigogine, 2015]). Brie y, we consider small deviations from this homogeneous solution (δc A ) and analyse them in Fourier space to nd if these perturbations decay back to the homogeneous solution or they are ampli ed and the system becomes unstable. To linear order in δc A and δc B , in Fourier space we have where k is the Fourier wave number and δµ

∂δc A ∂t = -k 2 D A δµ A + k 1 (e β(µ B +∆µ) δµ B -e βµ A δµ A ) + k 2 (e βµ B δµ B -e βµ A δµ A ) (5.6) ∂δc B ∂t = -k 2 D B δµ B -k 1 (e β(µ B +∆µ) δµ B -e βµ A δµ A ) -k 2 (e βµ B δµ B -e βµ A δµ A ), (
A = δc A 1 c A -2χ + κk 2 , δµ B = δc B 1 c B + κk 2
. Solving this linear system, amounts to nding the eigenvalues and eigenvectors of the Jacobian matrix; the stability of the homogeneous solution being ensured if, for every k, both eigenvalues are negative.

The sign of the determinant encodes information about the sign of the eigenvalues, since it is equal to their product. Thus, in order to assess the stability of the homogeneous system, we analyse the sign of the determinant of the Jacobian matrix (J) of the system de ned by Eqs. (5.6) and (5.7). To uncover the behaviour of the system, we simplify the super uous parameters as follows: D A = D B = 1. Then, for small wave numbers (k → 0, long-wavelength perturbation) the determinant takes the form

det(J) = 2(1 -2χc 0 A )(k 1 + k 2 e β∆µ ) c 0 A k 2 + O(k 4 ), (5.8)
where the superscript '0' indicates the value of the homogeneous concentrations. By inspection, one can see that the change of sign in the determinant corresponds to the spinodal condition for equilibrium systems (5.2), because k 1 and k 2 are strictly non-negative. When the determinant changes sign, one of the eigenvalues of J also changes sign, marking the threshold for stability, which matches the stability condition for the mixture at equilibrium. This is in line with the ndings of Chapter 4 where we saw that the Lyapunov functional of complex-balanced solutions is close to that of a system at thermodynamic equilibrium and cannot sustain di usion currents at steady state.

Non-complex-balanced behaviour

Finally, if we assume one of the reaction pathways considered above is sensitive to changes in the environment (as discussed in Section 3.2.3) complex balance is broken. We can describe this situation by the following version of the previous reaction-di usion equations

∂c A ∂t =D A ∇ 2 µ A + k 1 g(c A )[e β(µ B +∆µ) -e βµ A ] + k 2 (e βµ B -e βµ A ) (5.9) ∂c B ∂t =D B ∇ 2 µ B -k 1 g(c A )[e β(µ B +∆µ) -e βµ A ] -k 2 (e βµ B -e βµ A ), (5.10)
where we have included the function g(c A ) to represent the variation of the rate constant of the rst reaction with the environment of the reaction and, in particular, with the concentration of species A.

Carrying out a similar LSA with further simpli cations to neatly observe the e ect of breaking complex balance (k 1 = k 2 = 1), one nds, for the value of the determinant of the Jacobian matrix for small k det 5.11) where g (c 0 A ) is the derivative of g with respect to c A evaluated at c A = c 0 A . Thus, the sign of the determinant for small k is no longer set by the equilibrium spinodal decomposition as there is a contribution from the derivative of g, the function that describes how the rate of a chemical C 5. B

(J) = (1 + e β∆µ )(4 -8χc 0 A ) + c 0 A g (c 0 A )(1 -e β∆µ ) 2c 0 A k 2 + O(k 4 ), (
λ + c 0 A g ′ (c 0 A ) = 0 g ′ (c 0 A ) = 0.5 g ′ (c 0 A ) = 2.5 λ -
reaction changes with the environment. It should be noted that, in order to have the stability modi ed by non-equilibrium e ects as found here, two conditions are necessary: ∆µ = 0 and g (c 0 A ) = 0 (and not equal for every reaction in the network). The rst one implies that the system is not at thermodynamic equilibrium and the second one gives a local measure of how the rate of one of the chemical reactions changes due to the environment it take place in. If any of these conditions is not respected, then the stability of the system would be determined by the equilibrium spinodal condition, as for the complex-balanced or equilibrium case.

Note that here we have only analysed the stability of the system for small k or, in other words, long wavelength perturbations. In order to examine it with more generality, numerically, the value of the largest eigenvalue (with any wavelength) for each parameter set was obtained, see Fig. 5.1. As it can be seen from the gure, the complex-balanced (g (c 0 A ) = 0, solid blue line) eigenvalues yield the same condition for the instability as for spinodal decomposition at equilibrium λ -, Eq. (5.2). This suggest that more exotic behaviour could appear for non-complex-balanced solutions, such as Ostwald ripening arrest for non-complex-balanced systems, as found in Ref. [START_REF] Kirschbaum | Controlling biomolecular condensates via chemical reactions[END_REF].

Therefore, we have given minimum conditions for non-equilibrium patterning to appearcomplex balance must be broken -but not a su cient condition. Further e orts will be required in order to give general and su cient conditions for non-equilibrium patterning to appear and characterise the phenomenology that may arise.

Perspectives

It is therefore necessary to go beyond complex-balancing in order to be able to explain nonequilibrium phenomena such as patterning, which are ubiquitous in living beings.

A small step in this direction could be the analysis of networks that we will call weakly noncomplex-balanced. This is a network whose topology is not limited to a complex-balanced network (that is, its de ciency δ is not 0). However, if the reactions that break the complex-balance 5.2. P 97 constraint are slow when compared to the rest, then one would expect that the complex-balanced behaviour still dominates the dynamics but, arguably, with a small deviation due to these slow reactions. However, how large could this deviation be and what would its functional form be remains unclear.

Hence, we suggest a series expansion around the complex-balanced steady state to address this problem. Plenty of series expansions have been devised to analyse the steady-state behaviour of systems out of equilibrium (recently, e.g. [START_REF] Freitas | Linear response in large deviations theory: a method to compute non-equilibrium distributions[END_REF]). Nevertheless, most of this expansions take as a starting point the state of the system at thermodynamic equilibrium which means that, in principle, they cannot account for phenomena happening far away from equilibrium. In contrast, in this case, the expansion would not be based in the equilibrium steady state but in the complex-balanced one, which can be arbitrarily far away from equilibrium. Thus, we would be limited by the topology of the network and the timescales of the reactions rather than by the distance to thermodynamic equilibrium.

In practice, this would mean to obtain the steady-state distribution for the closest complexbalanced network and then add perturbatively the e ects of the slow reactions that break complex balance, where the expansion parameter would be the ratio between the reaction constant of the fast reactions and the slow ones. This could be done in an analogous way to other expansions (e.g. [START_REF] Proesmans | Linear stochastic thermodynamics for periodically driven systems[END_REF], Freitas et al., 2021]) but taking as an ansatz the complex-balanced steady state plus a correction, that is, an exponential with argument β

[F (n) + i μi n i + g(n)],
where g(n) quanti es the deviation from complex balance and could be obtained to di erent orders in the ratio of reaction constants.

Apart from increasing our understanding of general reaction networks, this expansion would also have consequences from the practical point of view. In real biological systems, it is seldom the case that the full reaction networks are known with great precision. Therefore, it might occur that a CRN which, to the best of our knowledge, is complex balanced, in reality it is not so, because there could be slow or rare reactions taking place that are hard to measure or identify. Then, this expansion would provide a means by which the results presented here could be robustly applied to complex biological settings, where deviations from complex-balancing, even if unidenti ed, may occur.

Finally, it would still be a challenge to nd Lyapunov functionals even for the weakly noncomplex-balanced case, if they exist [START_REF] Ge | Mesoscopic kinetic basis of macroscopic chemical thermodynamics: A mathematical theory[END_REF]. This is because arguments like the ones given in Section 4.5 are unlikely to hold for more complex steady states. In fact, any deviation from a Boltzmann distribution other than a linear one at the level of the exponential of particle numbers (which is the case of complex-balanced networks) would invalidate the proof of the Lyapunov functional for the di usion process, at least as performed here. Thus, in the future, new arguments to nd Lyapunov functionals will need to be devised, in order to obtain a more general understanding of non-equilibrium patterning events. But, for the moment, to the best of my knowledge, these more general Lyapunov functionals remain unknown.

C

The present work has been devoted to the study of theoretical biophysics and, in particular, patterning and self-organisation processes in non-equilibrium solutions. This topic has been analysed from two di erent perspectives, embodied by each of the two major parts of this manuscript.

In Part I, we theoretically addressed two examples of processes that organise the intracellular space in di erent ways: one that sets the size and number of large protein aggregates and another one which explains the appearance of a nucleoid in bacteria and its localisation within the cell.

To a great extent, the scienti c knowledge is built upon confrontation between hypotheses and facts, or, in this case, between theoretical predictions and quantitative measurements. Proceeding in accordance with these ideas, analysis of microscopy images were performed to confront the aggregation model. We found good agreement between the model and the experimental data and hypothesised that such a timescale separation could in uence many other intracellular phenomena that involve large agglomerates or structures, such as phase-separated droplets. Regarding the modelling of the bacterial nucleoid, our results were compared with the existing empirical evidence and additional measurements were suggested for several perturbations of the bacterial nucleoid, which should further validate the theory or force us to modify it, or even fully discard it.

In contrast, the objective of Part II was markedly di erent. In a purely theoretical endeavour, its aim was to understand the behaviour of reaction-di usion models in the most general way possible, both from the dynamical-system perspective (for example, what kind of behaviour can we expect from the reaction-di usion equations at steady state) and from the thermodynamical point of view (e.g. how does thermodynamics constrain the form of the reaction-di usion equations). We achieved a detailed characterisation of a large class of CRNs known as complex-balanced, for which we ruled out the presence of di usion currents at steady state (and, hence, non-equilibrium patterning cannot exist for these networks) and obtained a Lyapunov functional for these systems with a clear thermodynamic interpretation. In doing so, results that had been obtained for speci c networks (e.g. [Bazant, 2013, Kirschbaum andZwicker, 2021]) were generalised, but the resulting theory is not general enough to explain how non-equilibrium patterning arises in non-ideal solutions.

However, it remains unclear what are the repercussions of the results obtained in Part II for complex biological systems. Clearly, this is the price one has to pay for constructing such an abstract and general theory: it is still far away from making concrete falsi able hypothesis, especially in relation to biological systems. The applicability of these results is also limited by their heavy dependency on the topology of the CRN, which is well de ned from a mathematical point of view but in real living beings it might be hard to obtain in full detail. In particular, let us assume that, for a certain process in a given organism, molecular biologists have managed to discover the main players and their interactions. However, one cannot expect that this network is completely isolated from the rest of the organism nor that all chemical reactions have been discovered, no matter how weak or slow they are. Then, if the known network were complex-balanced but there were still some hidden nodes that break this complex-balancing, in view of Part II, what can we say about the system? Chapter 5 hints at how progress in this direction may be achieved, which may enable the robust application of these results in biological systems.

Nevertheless, the utility of this theory does not lie as much in its predictive power as in the C 101 order that it introduces in the plethora of results that have been obtained for particular models: this abstract approach enables the identi cation of a large class of networks that, even if driven out of equilibrium, still have many of the features associated to systems relaxing to thermodynamic equilibrium. From the theoretical point of view, it aids in relating network topologies to the phenomenology of non-ideal solutions out of thermodynamic equilibrium and has the potential to help in the construction of future models because it has increased, in an intuitive way, our understanding of the behaviour that these models may exhibit.

The original motivation for the work undertaken in Part II is the existence of e ective reactiondi usion models for non-ideal solutions that, while they are not thermodynamically consistent in the sense of Part II (chemical reactions and di usion do not stem from a single free energy), have achieved a remarkable accuracy in describing experimental observations (e.g. Ref. [START_REF] Zwicker | Centrosomes are autocatalytic droplets of pericentriolar material organized by centrioles[END_REF]). Even the model for the nucleoid presented in Chapter 2 lacks thermodynamic consistency according to Part II, yet it is a simple model that can recapitulate plenty of the experimentally observed phenomena and produce new hypotheses, which may or may not be con rmed. Both are simple but e ective models that could be thought of as coarse-grained models, where the microscopic details have been integrated out or absorbed into certain parameters, and that display this paradoxical behaviour under certain approximations (in the sense of an apparent lack of thermodynamic consistency in the dynamics). To the best of my knowledge, there is no systematic way to obtain these e ective models from a coarse-graining procedure that starts from microscopic and thermodynamically consistent models. This class of active matter models are not limited to reaction-di usion systems [START_REF] Li | Non-equilibrium phase separation with reactions: a canonical model and its behaviour[END_REF] but also encompasses systems with non-reciprocal interactions [START_REF] Saha | Scalar active mixtures: The nonreciprocal Cahn-Hilliard model[END_REF] or self-propelled particles [START_REF] Tailleur | Statistical mechanics of interacting run-andtumble bacteria[END_REF]. In some cases, attempting a coarse-graining procedure may not even make sense (e.g. for self-propelled particles it seems an overwhelming task) but, in others, I expect it to greatly increase our understanding of systems out of thermodynamic equilibrium, if such a mechanism could be devised. From our results, it seems clear that inherently active patterning phenomena as the one described in Ref. [START_REF] Li | Non-equilibrium phase separation with reactions: a canonical model and its behaviour[END_REF], Ref. [START_REF] Zwicker | Centrosomes are autocatalytic droplets of pericentriolar material organized by centrioles[END_REF], or Ref. [START_REF] Saha | Scalar active mixtures: The nonreciprocal Cahn-Hilliard model[END_REF] cannot arise from the coarse graining of a complex-balanced CRN1 , but, by no means this implies that these e ective models cannot be obtained starting from more general networks.

Each of the parts of this thesis represents two orthogonal approaches to theoretical biophysics. A modelling approach to theoretical biophysics utilises e ective representations of processes while retaining their physical interpretation and eases the obtention of clear falsi able predictions. This enables rapid rejection of the model if the hypotheses are proven wrong and advancement in the theoretical and empirical understanding of the system under consideration. Conversely, we have also undertaken a more abstract approach, aiming at deciphering the generalities of a certain class of reaction-di usion systems, which enables the extraction of general features of a wide variety of models in an intuitive way. However, as we have seen above, the drawback of this type of general intuitive understanding is the di culty in testing such a theory in real systems and the limitations when trying to explain complex emergent phenomena like non-equilibrium patterning. C Therefore, to make further progress in the understanding of the physics underpinning biological systems, both approaches will be needed: general results that set the main lines of the behaviour of non-equilibrium systems and particular models that yield concrete refutable hypotheses regarding the functioning of biological cells or organisms. Ideally, eventually these two approaches would be linked by a rigorous coarse-graining procedure that, in an understandable way, explains the functioning and self-organisation of the living from physical principles. Tackling the complexity of the living and comprehending the emergent phenomena in biological systems will require theoretical e orts at every scale, and a successful integration of theoretical predictions and experimental results that enables the empirical testing of general and abstract theories. does not yield the phase separation in the presence of reaction and out-of-equilibrium terms, while the addition of the third virial term provides such phase separation.

The second virial term consists of a summation of terms like the following

1 V N k dx k f (x i , x j ) = 1 V dr ij e -βV ij (r ij ) , (B.4)
where we have used the fact that the system is one dimensional and integrated out all spatial variables except one, r ij = x i -x j . The potential V ij depends on the relative position of species i and j considered, and here is taken to be a hard-sphere potential, which is non-zero only when the distance between two particles is smaller than the sum of their respective radii (R i + R j ). Mathematically, it takes the following form:

V (r 1 , r 2 ) = 0 if |r 1 -r 2 | ≥ R i + R j ∞ if |r 1 -r 2 | < R i + R j . (B.5)
For a given species i there are N i (N i -1)/2 combinations of the Mayer f function: as a result

N a (N a -1) V dr ij e -βV ij (r ij ) = - N a (N a -1) 2 V 4π 3 (2R a ) 3 . (B.6)
Conversely, as far as the hard-sphere interaction between species A and B is concerned, we have N A N B terms, yielding

N A N B V dr ij e -βV ij (r ij ) = - N A N B V 4π 3 (R A + R B ) 3 ≡ - N a N b V B (2) AB . (B.7)
In general, the term B (n) ab , which enters the n-th order of the virial expansion for species a and b, is known as virial coe cient. With this new notation, the partition function to order O(f

(x i , x j )) is Z = Z ideal    1 -   a=A,B N a (N a -1) 2 V B (2) aa + N A N B V B (2) AB      , (B.8)
where the partition function of an ideal gas of indistinguishable particles is

Z ideal = V N a=A,B N a !Λ 3Na a , (B.9)
and Λ a is the thermal de Broglie wavelength of species a.

The third virial coe cient is obtained by adding up the contributions of the three-particle interactions in Eq. (B.1), i.e., terms of the form

dx i dx j dx k f (x i , x j )f (x j , x k )f (x k , x i ). (B.10)
Proceeding along the same lines as for the rst virial term, we count the number of times when such a term appears in Z, which depends on the chemical species of the particles involved. For clusters of three particles all of the same chemical species we have the following contribution:

- a=A,B N a (N a -1)(N a -2) 6 160π 2 R 6 a 3 V 2 = - a=A,B N a (N a -1)(N a -2) 6 V 2 B (3) aaa , (B.11) B.2. A 119
where the summation over a involves the di erent chemical species A and B. If two particles belong to the same species and a third particle to a di erent species, we have

a =b - N a (N b -1)N b 2 16π 2 9 V 2 (8R 3 a R 3 b + 15R 2 a R 4 b + 6R a R 5 b + R 6 b ) = a =b - N a (N b -1)N b 2 V 2 B (3)
abb .

(B.12) These contributions are added up to the free energy of the system, to which higher order corrections can be introduced by computing the integrals involving products of 4 or more Mayer f functions.

Taking the logarithm of the partition function (B.8), once the third virial correction has been added, yields the free energy (2.4) in the main text.

B.2 Auxiliary entropy

In our analysis, we rst determine the steady state of the system in the absence of reaction and out-of-equilibrium terms, by minimising the total free energy (2.20). These steady-state pro les will then be used as initial conditions to integrate forward in time the reaction-di usion Eqs. (2.1) to (2.3), which include both reaction and out-of-equilibrium terms. At the free-energy minimum, the DNA concentration is nonzero in the nucleoid, while it vanishes outside the nucleoid. Given that these equilibrium pro les are entered as initial conditions in Eqs. (2.1) to (2.3), the vanishing concentration above causes numerical instabilities when these equations are numerically integrated forward in time, and can lead to negative concentrations in the out-of-equilibrium steady state [Shampine et al., 2005]. To overcome this issue, we included a small, additional entropic term in the free energy (2.20): .13) where c DNA (x) = N DNA /(2 ) is the average DNA concentration across the cell, and K aux and D aux are constants that we set to 0.2 N DNA and 10, respectively, in order to obtain a negligible auxiliary free energy in the nucleoid and a non-negligible one at the poles. The exponential term in Eq. (B.13) is such that, if the coordinate x lies in the nucleoid bulk, then c DNA (x) is of the same order of magnitude as c DNA (x) and thus, by choosing D aux su ciently large, the contribution to F aux vanishes. On the other hand, the exponential is approximately equal to one outside the nucleoid, where c DNA (x) c DNA (x) , and Eq. (B.13) reproduces a term that resembles the standard entropy of mixing of a polymer in the mean-eld approximation In what follows, we will add F aux to the total free energy Eq. (2.20), by setting

F aux = K aux k B T - dx e -Dauxc DNA (x)/ c DNA (x) c DNA (x) N DNA log [2 c DNA (x)] , ( B
F → F + F aux . (B.15)
As a result, the minimization of Eq. (2.20) tends to also maximize the entropy (B.14), thus spreading out a fraction of DNA segments outside the nucleoid bulk. This procedure alters only slightly the free-energy minimum: As shown in Fig. B.1, the auxiliary free energy does not vanish outside of the nucleoid, but it is orders of magnitude smaller than the original free energy within the nucleoid. Notwithstanding this, such a small free energy prevents numerical instabilities in the integration of the reaction-di usion equations.

B.3 E ect of third-order virial terms

Finally, we can evaluate the e ect of adding a third virial coe cient in our model. In Fig. B.2, we depict the equilibrium minimum of the free-energy and the out-of-equilibrium steady state, to both second and third order in the virial expansion. There is a marked di erence, particularly for the out-of-equilibrium steady state, where, without a third order term, the phase-separated nucleoid disappears due to the non-equilibrium synthesis of mRNA.

B.4 Estimate of the di usion coe cient of DNA

In order to assess the validity of the value for the di usion coe cient of DNA segments used in the model, we estimated the drag force exerted on the nucleoid by the viscous uid surrounding it (which in the cell can be identi ed with the cytoplasm).

We consider a cylindrical nucleoid with radius R -L/2, assuming that there is a gap of width L/2 (the persistence length of DNA segments) in the radial direction between the nucleoid and the plasma membrane. For a nucleoid moving as a whole with velocity v nucl , we assume the radial velocity pro le of the viscous uid surrounding the nucleoid that decreases linearly from velocity v nucl at the nucleoid to zero at the membrane. We thus have a velocity gradient in the radial from the nucleoid to the membrane. Then, the shear force F nucl exerted on the nucleoid is [Lautrup, 2011] F nucl = A nucl ν ∆v, (B.17) where ν and A nucl are the viscosity of the uid and the area of the nucleoid, respectively. For the reference cell of length 1.8 µm, the nucleoid will be around 1µm long and have a radius of R -L/2, which gives an area for the nucleoid A nucl = 1.9 µm 2 . We use the viscosity of water ν ∼ 1 mPa • s, although the viscosity of the cytoplasm might be larger. The ratio between the force F nucl and the velocity of the nucleoid is de ned as the drag coe cient [START_REF] Kubo | Statistical physics II: nonequilibrium statistical mechanics[END_REF]]

γ = F nucl v nucl = 2A nucl ν L 4.7k B T s/µm 2 . (B.18)
Given that the chromosome of the reference cell is composed of N DNA segments, if one naively ignores the hydrodynamic coupling between segments, one could infer that the drag coe cient per segment is γ s = γ/N DNA . Then, by the uctuation-dissipation relation, the di usion coe cient of each DNA segment would be D = k B T γ s 10 3 µm 2 /s, (B.19) which is signi cantly higher than the value used in the model. While this estimate is clearly an oversimpli cation, the overall low drag of the nucleoid suggests that processes where the di usion 122 A B. T B N coe cient of segments could matter, i.e., nucleoid centering, are not limited by viscous drag (and, hence, by the di usion coe cient), but by the osmotic-pressure di erence.

Our actual choice of the di usion coe cient for the DNA segments is signi cantly lower than the value (B.19) obtained with this estimate to be consistent with the measured di usion coe cient of mRNAs that have a similar linear dimension while agreeing with the result of the estimate (B.18), namely, that viscous drag should not limit the dynamics of the nucleoid. Both in the centering dynamics and in the expansion of the nucleoid after halting transcription, the timescale is set by the synthesis and degradation of mRNA (β -1 ), and not by the di usion coe cient.

B.5 Experimental Methods

I include this section of the Appendix for completeness yet it is not my own work but that of my collaborators, as the experiments were carried out exclusively by Sophia Hsin-Jung Li.

In this study, we used E. coli wild type strain NCM3722. To achieve di erent growth rates, cells were cultured at 37 • C in chemostats and in batch. For slow growth rates (0.1 and 0.6 h -1 ), carbon-limiting chemostats with corresponding dilution rates were used, whereas for faster growth, batch cultures with glucose minimal media (0.9 h -1 ) and de ned rich media (1.7 h -1 ) were employed. The chemostat (Sixfors, HT) volume was 300 mL with oxygen and pH probes to monitor the culture. The pH was maintained at 7.2 ± 0.1 and the aeration rate was set at 4.5 l h -1 . 40 mM MOPS media (M2120, Teknova) was used with glucose (0.4 %, Sigma G8270), ammonia (9.5 mM NH4Cl, Sigma A9434), and phosphate (1.32 mM K2HPO4, Sigma P3786) added separately. For the de ned rich media, additional Supplement EZ 5X and 10X ACGU Solution (Teknova) were added. In carbon-limiting chemostats, glucose concentration was reduced to 0.08 %. All the sample collection happened after chemostat cultures reached steady state or when batch culture reached OD600 0.3.

To measure cell size, 750 µL of culture was xed with 250 µL 20% paraformaldehyde at room temperature for 15 min, washed with PBS twice, and stored at 4 • C until imaging. Then, 1 µL of cells were placed on 1% low-melting agar pad (Calbiochem) made with PBS and imaged with inverted Nikon90i epi uorescent microscope equipped with a 100 × 1.4 NA objective (Nikon) and Hamamazu Orca R2 CCD camera. NIS Elements software (Nikon) was used to automate image acquisition for phase contrast images. Segmentation, quanti cation of uorescence intensity, and cell-length measurements were further analyzed in MATLAB [MATLAB, 2018] using customized programs.

To infer ribosome number per cell, cell number per OD600 and total RNA per OD600 were measured separately. Cell number per OD600 was calculated by serial dilution and plating. To measure total RNA, 1.5 mL of culture was pelleted by centrifugation for 1 min at 1.3 × 10 4 X g. The pellet was frozen on dry ice and the supernatant was taken to measure absorbance at 600 nm for cell loss. The pellet was then washed twice with 0.6 M HClO4 and digested with 0.3 M KOH for 1 hour at 37 • C. The solution was then precipitated with 3 M HClO4 and the supernatant was collected. The pellet was re-extracted again with 0.5 M HClO4. The supernatants were combined and absorbance measured at 260 nm using Tecan In nite 200 Pro (Tecan Trading AG, Switzerland). Total RNA concentration was determined by multiplying the A260 absorbance with 31 B.6. N 123 (µg RNA/mL) as the extinction coe cient.

B.6 Numerical methods

For both the minimisation of the free energy (2.20) and the time integration of the reactiondi usion Eqs. (2.1) to (2.3), numerical methods are required, as the equations are too complex to be solved analytically. We discretise the spatial degrees of freedom of the system into a mesh that satis es: (B.20) where ∆x is the distance between two neighboring points and N d is the number of points taken to describe the concentration pro le of each of the species in the system, which was set to N d = 32, 64, depending on the desired accuracy.

2 = ∆xN d ,
Using this discretisation of space, we evaluated the spatial derivatives in Eq. (2.20), and obtained a minimum of the free energy by using an algorithm for constrained gradient-based optimization [Kraft, 1994]. We used the C implementation of the NLopt library [Johnson, 2007].

For the time integration with the non-equilibrium e ects, we wrote down a set of ordinary di erential equations, where to each chemical species and each point in the mesh (de ned by Eq. B.20) corresponds a function of time, and such functions are coupled to neighboring points in space through the discretized spatial derivatives, and to other chemical species through the local chemical reactions. To solve this system we used an implementation of the backward di erentiation formula (BDF) method in Mathematica [Wolfram Research, Inc., 2018].

B.7 Parameter estimation from experimental data

Cell length, cross-sectional radius and number of ribosomes were inferred from experimental measurements of E. coli colonies growing in di erent chemostatted conditions (as explained in Section B.5). The data yields the values of these parameters for di erent growth rates. The estimate of the number of ribosomes is derived from the total amount of 23S and 16S ribosomal RNA (rRNA), considering that two-thirds of the total mass of a ribosome comes from rRNA [START_REF] Berg | Biochemistry[END_REF].

We analyzed twelve di erent nutrient limitations, which correspond to four groups of similar growth rates: the parameters-cell length, cross-sectional radius, number of ribosomes, and growth rate-were averaged across data points belonging to the same growth rate.

Finally, the values of the parameters above as functions of the growth rate were tted with an exponential by using the least-square method. The parameter values for a growth rate of log(2)/2/hr, which corresponds to a doubling time of 2 hr, were obtained via interpolation, by evaluating the exponential t at the reference growth rate log(2)/2 hr -1 , see Section 2.

B.8 Nucleoid size at thermodynamic equilibrium

In this section we depict the results obtained for the lamentous growth scaling as in the main text, Fig. 2.4, but in the absence of the out-of-equilibrium chemical reactions, i.e., for a passive system (see Fig.

B.9 Single chromosome growth condition B.9.1 Scaling of the concentration of chemical species for single-chromosome growth

In the single-chromosome case we assume that the concentration of mRNA and ribosomes scales linearly with the growth rate and, based on the data of Ref. [Kohram, 2021], that the growth rate of E.coli decreases linearly with cell length until it reaches zero at ∼ 20 µm. In order to compare the model predictions with the experimental data in Ref. [Wu et al., 2019], we assume the same linear law, but with a slope such that the growth rate, g, reaches zero at 30 µm, because in the data of Ref. [Wu et al., 2019] cells appear to grow up to that length. In addition, we assume that the mRNA degradation rate, β, decreases linearly in the same way the growth rate does, motivated by the expectation that for slow growth it would be ine cient to turn over mRNA quickly.

With the above considerations, we can write the following relations:

B.9.2 Results for constant RNA concentration

This gure depicts the predicted nucleoid length for single-chromosome growth, as in Section 2.3.2 Single-chromosome growth, but keeping mRNA and ribosome concentration constant. As it can be seen in Fig. B.5, the predicted nucleoid length is much smaller than the measured one, which led us to conclude that mRNA and ribosome concentrations will fall with cell length in single-chromosome lamentously growing cells, which is consistent with the observed reduction in growth rate in Ref. [Kohram, 2021].

A C. N -CRN di erence of the complex:

F (n) -F (n -r ρ ) = ∆F id + α r ρ α β χ α,β V n β - α,β χ α,β 2V r ρ β r ρ α + α χ α,α 2V r ρ α (C.4)
where ∆F id is the ideal part of the free energy di erence of the complex, as given by eq. (3.14). In the last line of eq. (C.4), only the rst term in χ α,β scales with the particle number, hence, in the thermodynamic limit, the rest of the interacting terms are negligible. In general, the propensity function in the stochastic system for an interacting regular solution is complicated, as seen in Eq. (C.4), but for a unimolecular reaction it simpli es to:

f +ρ (n) = k ρ e β(F (n)-F (n-r ρ )) = k ρ e βµ i = k ρ n i V e β(µ 0 i + j χ i,j V n β ) , (C.5)
where i is the reactant of the reaction +ρ and the de nition of chemical potential (3.18) has been used.

In the thermodynamic limit, the deterministic rates can be written, according to eqs. (3.12) and (C.4) as follows We note that reaction rates of this form have already been suggested in this deterministic limit in the context of phase-separating battery materials [Bazant, 2013] and that a similar result for the reaction rates can be derived from activity coe cients (see Appendix C.2), which are common in the chemistry literature.

J +ρ = k +ρ α c r ρ α α exp   α,β

C.2 Form of the reaction rates in the non-ideal solution theory

Chemists have been dealing for long with non-ideal solutions, such as the one we are attempting at describing, and its chemical potential its typically written as [Kondepundi and Prigogine, 2015] µ α = µ 0 α + log(γ α (c)c α ), (C.9)

where γ α (c) is the activity coe cient which, itself, can depend on the vector of concentrations of the mixture c. Note that, if we identify γ α (c) = exp β χ α,β c β we recover the regular solution chemical potential in the thermodynamic limit eq. (C.8). Then, the reaction rate becomes where we have made explicit the height of the barrier, µ * , and implicit the standard chemical potentials of the reactants. µ * can be interpreted as the free energy of the activated complex (cf. Fig. 3.1) and can be identi ed with µ * 0 + log γ * , that, after absorbing µ * 0 into the reaction constant yields the rate [START_REF] Madon | Catalytic reaction rates in thermodynamically non-ideal systems[END_REF]] (C.11) This expression of the reaction rates from the point of view of non-ideal solutions has the advantage of providing an intuitive explanation for the e ect of changes in concentrations on the reaction rates. Changes in concentrations will a ect the rates not only through changes in the chemical potential of reactants and products but also, according to eq. (C.11), through the degree of solvation of the activated complex * [START_REF] Madon | Catalytic reaction rates in thermodynamically non-ideal systems[END_REF], parametrised by γ * . This interpretation of the reaction rate gives an expression with a physical meaning for how the reaction rate should vary with the environment in which the reaction takes place.

R +ρ = k+ρ γ * α (c α γ α ) r ρ α .

C.3 Minimisation of the Lagrangian to obtain the phase diagram

In order to nd the steady state of the solution we need to minimise the Lyapunov functional (or the Lagrangian, once particle conservation constraints are taken into account), Eq. (4.47). Substantial simpli cation is made by neglecting the contribution of the interfaces and considering the system as two homogeneous phases, which implies that the actual function that needs to be minimised is the Lagrangian (4.48).

First, we reduce the dimensionality of the problem by equating the derivatives of the Lagrangian with respect to the concentrations of the species: where i and i are any two given chemical species in the system, and (k) refers to the di erent phases in the system. In practice, this yields relationships of the analogous to the equality of chemical potentials at equilibrium (which here include the shifted chemical potential term μ due to the non-equilibrium complex-balancing). For simple free energies like Eq. (4.49), we have

φ B =φ A e 2χφ A +μ A -μ B , (C.13) φ C =φ A e 2χφ A +μ A -μ C , (C.14)
which reduces the problem to just three variables: φ

(1)

A , φ

(2)

A and V (1) . We nally need to enforce stationarity condition of the Lagrangian with respect to volume: A , φ

(2)

A and V (1) . However, the resulting conditions are transcendental equations which, in general, have no explicit analytical solution. Therefore, they need to be solved numerically. Even numerically, it is still a hard problem for parameter sets near the critical point, which is why in Fig. 4.6 the density of data-points around the critical point decreases.

MOTS CLÉS

Physique statistique hors-équilibre, biologie cellulaire, organisation spatiale, démixtion, solutions non-idéales, réseaux de réaction chimiques.

RÉSUMÉ

Le cytoplasme d'une cellule est constitué d'une myriade de macromolécules comme des protéines des lipides et des enzymes. Chacune de ces molécules jouent un rôle bien spécifique au sein de la cellule, ce qui requiert un contrôle précis de leurs position spatiale. Dans la première partie de cette thèse nous proposons une explication théorique à deux processus où le milieu intracellulaire est organisé spatialement sous forme de motifs. Ces processus découlent tous deux du comportement collectif de ces molécules. Nous développons dans un premier temps un modèle qui a pour but d'étudier l'impact de l'encombrement stérique intracellulaire sur la cinétique d'agrégation de ces protéines. Nous montrons ensuite que les prédictions de ce modèle corrèlent de manière satisfaisant avec les données expérimentales de nos collaborateurs. Dans un deuxième temps, nous étudions la compaction et la localisation spatiale du chromosome de la bactérie Escherichia coli (E. coli). Nous expliquons la condensation du chromosome par les répulsion stériques entre l'ADN et d'autres macromolécules comme les ARNs messagers ou les ribosomes. Cette transition de phase s'assimile d'un point de vue théorique à la démixtion entre deux liquides. La localisation particulière du genome d'E. coli est quand-à-elle expliquée par l'activité transcriptomique des ARNs messagers. Cette hypothèse nous permet de reproduire la signature spécifique du chromosome d'E. coli in vivo, situé à la moitié de la cellule avant la division du chromosome et au 1/4 et au 3/4 après sa division. Dans la deuxième partie de cette thèse, nous construisons un cadre thermodynamique pour décrire les réseaux de réactions chimiques au sein de solutions non idéales. Cette approche nous permet de généraliser certains résultats de la théorie des réseaux de réactions chimiques idéales. De plus, ce cadre aide à formaliser les connections entre les réactions chimiques hors-équilibre et les séparations de phase pour les réseaux de type «complex-balanced». Ces derniers offrent un cadre moins restrictif à la théorie des réseaux de type «detailed-balanced». Nous caractérisons les réseaux de type «complex-balanced» par leur topologie. Nous montrons alors que cette topologie peut contraindre la dynamique des solutions. À l'inverse, notre théorie prédit que lorsque le réseau n'est plus de type «complex-balanced» des dynamiques plus exotiques peuvent apparaître.

ABSTRACT

In the cytoplasm of a biological cell there are a myriad of different proteins, lipids and enzymes, each of them performing different tasks. The spatial organisation of these chemical species is crucial for the correct functioning of a cell. In the first part of this thesis we will explore, from a theoretical perspective, two processes where the intracellular medium is patterned and organised, whose common feature is the fact that they both stem from the collective behaviour of a large number of molecules. First, we develop a model for protein aggregation which studies the effect of intracellular obstacles on the coagulation kinetics. Our predictions are then successfully compared with experimental data obtained by our collaborators. The other example refers to the compaction and location of the bacterial chromosome. We suggest that the chromosome segregates from the rest of the cytoplasm, because of steric interactions between DNA and the intracellular crowders, by means of a mechanism reminiscent of liquid-liquid phase separation. Moreover, our study indicates that spatial localisation of the chromosome within the cell is dictated by non-equilibrium transcription of mRNAs (which are part of the crowding effect). Our model successfully reproduces the localisation pattern of the Escherichia coli chromosome, which is positioned at the center of the cell before the division of the chromosome and at 1/4 and 3/4 of the cell after division. Building on these examples, in the second part of the thesis we construct a thermodynamically consistent framework to mathematically describe chemical reaction networks in non-ideal solutions. This framework allows us to generalise the results from the classical theory of ideal networks, and aids in elucidating the connection between non-equilibrium chemical reactions and phase separation for a large class of networks, known as complexbalanced networks. Complex-balanced networks are a class of chemical reaction networks that is less restrictive than detailed-balanced ones. Given that complex balancing can be fully determined from the topology of the network, we discuss how this topological property of the network can constrain the dynamics of the solution and what behaviour one can expect when complex balancing is broken.

KEYWORDS

Non-equilibrium statistical physics, cell biology, spatial organisation, phase separation, non-ideal solutions, chemical reaction networks.

  size for single-chromosome lamentous cells with constant mRNA concentration . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 126

Figure 1 . 1 :

 11 Figure 1.1: Experimental system studied and theoretical model. (A) Fluorescent images of a cell expressing CRY2olig and activated with blue light every 2 min. The cell shows many small clusters by t = 3 min, which then mature and coalesce over time. (B) Cartoon depicting the main ingredients of the model. Left: a density ρ of monomers that can freely di use and aggregate fast. Middle: N (t 0 ) larger clusters of size ∼ m * are trapped in the cytoskeleton. Right: N (t f ) large clusters are only able to move and further aggregate when the remodelling of the cytoskeleton sets them free.

  Figure1.2: Experimental data and ts to the theoretical model. Results obtained for di erent cells 1 hr after the beginning of the aggregation process. In (A) and(B) we plot the cluster concentration and mean cluster size as functions of the protein density ρ measured in each cell, respectively. Grey dots correspond to results for individual cells, black squares to the average over 20 cells, and error bars to standard deviations. Blue lines correspond to least-square ts of the theoretical expressions for the cluster concentration and cluster size, that is, Eqs. (1.28) and (1.29). The orange dashed lines correspond to the predictions for a passive cytoskeleton, i.e., in the absence of an active dynamics (β = 0). (C) Cluster-size probability density function (PDF) for di erent protein concentrations: Low (below 300 monomers/µm 3 ), medium (between 300 and 600 monomers/µm 3 ) and high (above 600 monomers/µm 3 ). The black dashed line corresponds to the estimate of the cluster-mass threshold between timescales, m * , inferred from (A) and(B).

Figure 2 . 1 :

 21 Figure 2.1: Sketch of the model for the nucleoid. (A) Cartoon of an E. coli cell and its transcriptional-translational machinery, where the horizontal axis is the single dimension we consider. Blue coils represent mRNAs in polysomes, the red coil denotes the DNA plectoneme and ribosomes are shown in black. Blow-ups: polysome composed of an mRNA and n ribosomes with gyration radius R n (left) and plectonemic structure of the DNA with persistence length L (right). (B) Schematic of the components underlying the reaction-di usion model, where the DNA plectoneme is represented by a set of disjoint cylinders, and polysomes and free ribosomes by spheres.

  Figure 2.2: Steric interactions of DNA segments.Hard-core interactions between DNA segments (red cylinders) and ribosomes or polysomes (black spheres). DNA segments interact with each other through a cylinder of radius ρ , while they interact with polysomes and ribosomes with a radius ρ > ρ , satisfying Eq. (2.13).

  a and b. It follows that κ ab is an intrinsic feature of the particles of species a and b:

  Fig. B.1.

  ρ ′ ) 2 and DNA segments [c DNA,tot , to distinguish from the position dependent DNA density c DNA (x)],

n

  Nucleoid length (µm) Cell length 2ℓ (µm) ρ tot = 1500 /µm Nucleoid length (µm) Cell length 2ℓ (µm) Nucleoid length (µm) Cell length 2ℓ (µm) ρ tot = 2400 /µm Nucleoid length (µm) Cell length 2ℓ (µm) Nucleoid length (µm) Cell length 2ℓ (µm) ρ tot = 3000 /µm Nucleoid length (µm) Cell length 2ℓ (µm) Nucleoid length (µm) Cell length 2ℓ (µm) Nucleoid length (µm) Cell length 2ℓ (µm) Nucleoid length (µm)

Figure 2 . 4 :Figure 2

 242 Figure 2.4: Steady-state concentration pro les for E. coli growing lamentously. (A) Concentrations along the long axis of the cell of DNA, c DNA (x) (red), and free ribosomes, c F (x) (black), and polysomes ρ n (x) with the mRNA loading number n indicated by the colour bar.Each panel corresponds to a di erent cell half-length , marked on the top of the panel, with total mRNA density ρ tot = 2400µm -1 . (B) Nucleoid length versus cell length, 2 , for di erent values of the total mRNA density, ρ tot . For each value of ρ tot the nucleoid length is shown up to the cell length at which the nucleoid splits into two lobes. Red data points were obtained from the numerical solution of the entire model while the grey lines are the analytical estimates. The inset depicts on a log-log plot the power-law relationship between the length at which the nucleoid splits and the rate of mRNA synthesis αc DNA , where c DNA is the average DNA concentration along the nucleoid at steady state. The triangle represents the scaling obtained from the simpli ed model of Eq. (2.49).

Figure 2

 2 Figure 2.6: mRNA uxes at steady-state.The plot depicts the concentration proles at the non-equilibrium steady state. The blue arrows in the bottom depict the steady-state ux of mRNAs, J ρtot (the arrow length is proportional to the mRNA current). At steady state, mRNAs are synthesised within the nucleoid, wherein they di use before escaping to the poles (and are ultimately degraded).

Figure 2 . 7 :

 27 Figure 2.7: Transcription splits the nucleoid into two lobes located at 1/4 and 3/4 of the cell length. (A) Concentration pro les for a lamentous cell, obtained from the equilibrium pro le at t = 0 by integrating forward in time the reactiondi usion Eqs. (2.1) to (2.3) in the presence of the non-equilibrium processes until t = 20 min and t = 40 min, for a cell with a half-length = 4.95 µm shown as in Fig. 2.4 A. (B) Positions of the center of mass of the left (dashed red curve) and right (solid red curve) halves of the DNA along the long axis of the cell as a fraction of the total cell length, as functions of time.

  Figure 2.8: Non-equilibrium processes centre the nucleoid at midcell. (A) Concentration pro les obtained by initially shifting the steady-state pro les towards the right cell pole at t = 0, and then integrating forward in time the reaction-di usion Eqs. (2.1) to (2.3) in the presence of nonequilibrium processes to t = 5 min and t = 20 min, for a cell with half-length = 1.8 µm, shown as in Fig. 2.4 A. (B) In red, location of the center of mass of the nucleoid along the long cell axis, as a function of time. In gray, the analytical lower bound obtained by neglecting nucleoid drag. Inset:The quantities depicted are the same as in B, but with the y-axis is in logarithmic scale.

Figure 2 . 9 :

 29 Figure 2.9: The nucleoid expands in the absence of mRNA synthesis. (A) Steadystate pro le including mRNA synthesis (t = 0) and pro les obtained by integrating forward in time from the steady-state pro le at t = 0 in the absence of mRNA synthesis (t = 4 mins and t = 20 mins), for a cell with half length = 3.6 µm. The concentration pro les are shown as in Fig. 2.4 A. (B) Fraction of the cell volume occupied by the nucleoid, as a function of time (computed as the fraction of length along the axis of the cell with a DNA segment concentration c DNA > 1000 µm -1). The turquoise data points were obtained from Cabrera et al.[START_REF] Cabrera | Active transcription of rRNA operons condenses the nucleoid in Escherichia coli: Examining the e ect of transcription on nucleoid structure in the absence of transertion[END_REF], for the case of a cell treated with rifampicin (which blocks mRNA transcription).

Figure 2 . 10 :

 210 Figure 2.10: The nucleoid contracts in the absence of translation initiation. (A) Steady-state pro le including translation initiation (t = 0) and pro les obtained by integrating forward in time from the steadystate pro le at t = 0 in the absence of ribosome binding (k on = 0 for panels t = 2 mins and t = 10 mins). (B) Time evolution of the fraction of the cell volume occupied by the nucleoid, shown as in Fig. 2.9B.

Figure 2

 2 Figure2.11: Steady-state pro les for lamentously growing cells. Steady-state pro les of the concentrations of the components of the E. coli transcriptional-translational machinery (TTM) for lamentously growing cells with mRNA concentration ρ tot = 2400 µm -1 . Colours represent di erent chemical species of the TTM, as in Fig.2.4 A, and each of the panels corresponds to a di erent cell length (as indicated within each panel). To enable comparison between di erent panels, the lengthscale in all panels was kept constant. As in Fig.2.7, for cells larger than ∼ 8 µm the nucleoid splits in two lobes. Furthermore, for cells around ∼ 16 µm the nucleoid has three distinct lobes, suggesting that this may be a pattern with a characteristic length that exists also for longer cells.

Figure 3 . 1 :

 31 Figure 3.1: Reaction free energy landscape for a chemical reaction in two di erent phases. The horizontal dimension is the reaction coordinate and the vertical one speci es the height of the reaction free energy landscape F 0 (which is equal to F except for the fact that it does not include the entropic term

  Fig 4.1).

Figure 4 . 1 :

 41 Figure 4.1: Graphical representation of CRNs discussed in Section A hierarchy of steady states. The CRN has 6 complexes each of them represented in one of the nodes of the network: A, B, C, A + D, E and B + D. The ve reactions present in the CRN are numbered.

  Figure 4.2: Graphical representation of CRNs discussed in Example 3.The left CRN has de ciency δ = 1 and the right one is zero-de cient δ = 0. Therefore, the network on the right will necessarily have a complex-balanced steady state but this is not the case for the left-hand one.

Figure 4 . 3 :

 43 Figure 4.3: Graphical representation of the CRN used in Examples 4and 5. These two CRNs are necessarily complex balanced (δ = 0 in both cases). Unimolecular networks (as the one on the right) are always complex balanced. In the case of the left-hand side network, despite the presence of a bimolecular complex A+D, the network is still complex balanced, since all the A or D particles are created or destroyed through the complex A + D, hence, in this case, the steady-state condition implies the complex-balance condition (4.9). The left CRN corresponds to the one used in Example 4; the right one corresponds to Example 5.

Figure 4 . 4 :

 44 Figure 4.4: Steady-state probability distributions for Example 3. Marginal probability distribution (A)and conditional probability distribution(B) for Example 3 at steady state. The blue triangles are the probabilities obtained from a Gillespie simulation and the red line is the analytical solution given by Eq. (4.26).

ρi

  and b = i (µi + μi)s ρ i (and conversely for the second).

Figure 4 . 5 :

 45 Figure 4.5: The Lyapunov function of a complex-balanced CRN is minimised by its dynamics. Numerical results for the CRN (4.44) obtained from perturbing a homogeneous solution and integrating Eq.(4.38) forward in time until it reaches a steady state. Top panel: evolution of the value of the Lyapunov functional with time. This functional is minimised with the dynamics but it does not attain the value 0 for two reasons: First, the constant termlog(Z) has been neglected. Second, the dynamics may not have reached the global minimum of L but only a local one. Middle panel: Steady-state con guration of the system, displaying the volume fractions of each species as a function of the spatial coordinate φ(x). In this case, the system reaches at a steady state where there is phase separation. Lower panel: Net reaction ux at steady state for the third reaction in the CRN (4.44) as a function of the spatial coordinate, showing that the system is not detailed-balanced. Note that, while in the upper panel, the horizontal dimension refers to time, it is no longer the case for the lower ones, where it represents the space in a one-dimensional system.

Figure 4 . 6 :

 46 Figure 4.6: Phase diagrams for a chemically reactive non-ideal mixture. Lower panels: Phase diagrams obtained by minimisation of Eq. (4.48) as a function of total solute fraction φ N and interaction parameter χ or non-equilibrium driving ∆µ (as de ned in Example 5). The color code indicates the amount of A particles along the phase coexistence lines (volume fraction). Upper panels: Each of the panels has associated a number which states their corresponding point in the phase diagrams. These panels depict the steadystate volume fraction pro les of the system in a one dimensional space, x, from 0 to L. These results were obtained by numerically integrating in time the reaction-di usion equations, as in Fig.4.5, with the parameters set by their point in the phase diagram. Note that, in the upper panels, the width of the interface is nite to ease the numerical integration, unlike in the lower panels, where the phase diagram was obtained for vanishing interfaces.

Figure 5 . 1 :

 51 Figure 5.1: LSA for complex-balanced and non-complex-balanced networks.Largest eigenvalue (evaluated at the wave number k that maximises it) for each of the values of a homogeneous c 0 A . The eigenvalues λ + are compared with the equilibrium spinodal condition λ -(solid green line). λ -crosses 0 at the same point where the g (c a ) = 0 line becomes positive (solid blue line), proving that for complexbalanced systems, despite being out of equilibrium, their stability is still given by the equilibrium spinodal condition. The parameters are set as in Section 5.1, with ∆µ = 5k B T , c 0 B = 1 and g(c 0 A ) = 1.

c

  DNA (x) N DNA log [2 c DNA (x)] , (B.14)which dominates the integral in F aux .

  Figure B.1: DNA and auxiliary free energies. DNA and auxiliary free energies (F DNA and F aux , respectively)for a cell of size 3.6 µm at the out-of-equilibrium steady state. The DNA free energy only takes into account the self-interaction terms of DNA cylinders, because the mean-eld entropic term of the DNA chain is neglected, see Section 2.2.2.

B

  Figure B.2: E ect of third-order virial term. (A) Pro les corresponding to the free-energy minimization (top) and out-of-equilibrium steady state (bottom), obtained with the second virial coe cients only, for a cell length of 3.6 µm for single-chromosome growth. As in the main text, concentrations of DNA, c DNA (x), and free ribosomes, c F (x) are shown in red and black respectively and polysome concentrations ρ n (x) are also shown, where the mRNA loading number n is speci ed by the color box. (B) Same quantities as in A, with the third virial coe cients.

  Figure B.3: Experimental data and interpolation.(A) Ribosome number per cell as a function of growth rate. Points correspond to experimental data, and solid curves to the exponential t. The inferred parameters for the desired growth rate (log(2)/2hr -1 ) are marked on each axis. (B) Same as A, for cell cross-sectional radius. (C) Same as A, for half cell length.

  B.4). The format of the Figure is the same as that of Fig.2.4.

Figure B. 4 :

 4 Figure B.4: Equilibrium concentration pro les for E. coli growing lamentously. This gure shows the concentration pro les and nucleoid lengths exactly as in Fig.2.4 except that the non-equilibrium terms in the equations have not been added. Therefore, this gure represents the results of the model at thermodynamnic equilibrium.

  particle numbers have been replaced with continuous concentrations and the internal energy part of the chemical potential has been absorbed by the rate constant k +ρ . However, the rates (C.6) can be written in the more general formJ +ρ = k ρ exp α

R

  +ρ = k +ρ α (c α γ α ) r ρ α e -βµ * ,

The work presented in this chapter has been published in Biophysical Journal, Ref.[Miangolarra et al.,

2021a]. The present chapter is therefore largely based on that reference. The experiments described in this chapter were performed by Aléria Duperray-Susini and Mathieu Coppey.

The work presented in this chapter has been published in Proceedings of the National Academy of Sciences of the USA, Ref.[Miangolarra et al.,

2021b]. The present chapter is therefore largely based on that reference. The work was done in collaboration with Ned S. Wingreen and Sophia Hsin-Jung Li.21

In the literature, there is a large variability in the value for the persistence length of the DNA plectoneme, from the

30-50 nm suggested by[START_REF] Verma | Architecture of the Escherichia coli nucleoid[END_REF], Xiang et al., 2021] to the 75-100 nm used by[Odijk, 2000, Mondal et al., 2011]. We take 100 nm, because this is the value used in the simulations of Ref.[START_REF] Mondal | Entropybased mechanism of ribosome-nucleoid segregation in E. coli cells[END_REF].

From the previous section, e.g. Fig 2.4, one can check that this is a good approximation. From the free-energetic point of view, it makes sense to have a compartment enriched in DNA and depleted in crowders, since it decreases the interaction energy due to the form of the inter-particle potentials considered.

Even a completely random positioning within a cell would yield midcell as the average position.

This distinction has been made to avoid confusion with the jargon used in the mathematical literature where an equilibrium point in a dynamical system refers to a point in the system with vanishing time derivatives, independently of any thermodynamic considerations.

The diluteness condition has been introduced to ensure that, for the moment, interactions between solutes will not be strong enough to cause inhomogeneities in the solution. Afterwards, we will see that phase separation can occur in the presence of interactions between solutes and, thus, the system cannot be described by a single homogeneous concentration.

We have used a canonical Boltzmann distribution for the argument, which implicitly assumes that total number of particles xed. However, a generalisation can be obtained for a system in equilibrium with a particle reservoir by using the grand-canonical distribution instead.

Note that this is not the most general form of the propensity functions. A more general expression will be presented later in Section 3.2.3.

This expression di ers from the one given in Ref.[Anderson et al.,

2010] due to the explicit consideration that every reaction can happen in both senses, to better illustrate the alignment with thermodynamic principles.

Note that, since patterning can also emerge at thermodynamic equilibrium, we are only referring to patterning that arises as a consequence free-energy dissipation and that it would not exist in the absence of this dissipation.

Acknowledgements

Appendix A

Supplementary Information for the aggregation model A.1 Experimental Procedures

I include this section of the Appendix for completeness yet it is not my own work but that of my collaborators, as the experiments were carried out exclusively by Aléria Duperray-Susini and Mathieu Coppey.

A.1.1 Cell culture

The immortalized hTERT RPE1 cells (Human Retinal Pigmented Epithelium) were cultivated in DMEM F12 without Phenol Red (Gibco, Life Technologies) supplemented with 10% Fetal Bovine Serum (FBS) without antibiotic, hereafter called the growth medium. They were maintained at 37 • C in humidi ed atmosphere with 5% CO 2 , tested and certi ed as mycoplasma free.

A.1.2 Transitory cells transfection by Cry2Olig-mCherry

RPE1 cells were detached by trypsin and centrifuged for 3 min, 100 g at room temperature to eliminate it. The pellet was kept and resuspended on growth medium. They were transfected on suspension by jetPrime (Polyplus transfection), with 1 µg of DNA plasmid vector Cry2Olig-mCherry (purchased from Addgene, number 60032), and then platted on uorodishes. According to the recommendation of manufactory, the medium was replaced after four hours by a fresh one. From there, the manipulation of cells was done in the complete dark.

A.1.3 Quantitative estimation of uorescent protein concentration

To estimate the concentration of proteins in cells using the uorescent signal, we calibrated the intensity on the camera using mCh-6His protein puri ed at 4.19 mg/ml (a gift from El Marjou. A, Platform of Curie Institute). We performed serial dilutions of the stock solution (1, 1:2, 1:4, 1:8, 1:10, 1:16, 1:32, 1:64, 1:100, 1:128, 1:1000) in the cell growth medium, and the medium alone was used for background estimation. For each dilution, we put a drop of 10 µl into a uorodish 113 A A. A P and we imaged the drop using the exact same parameters as for the cell imaging experiments. Two images were acquired at a focus right above the coverslip, as for cell imaging. We then quanti ed the average uorescent intensity using Fiji. The total intensity of the image was background subtracted and averaged over the size of the whole image. Data were plotted and gave rise to a linear relationship between raw intensities of the images and concentrations of recombinant uorescent proteins. We tted data with a line and used the value of the slope to convert intensities into concentrations.

A.1.4 Optogenetic experiments

All experiments were performed using 100x objectives (oil immersion, numerical aperture 1.4) by Inverted Spinning Disk Confocal Roper/Nikon, EMCCD 512x512 evolve (pixel size: 16 µm) photometrics come from to Imaging Nikon Center (PICT-LM) in Curie Institute. Live imaging was on normal growth condition and preserved by Life Imaging Service Yokogawa head: CSU-X1 integrated in Metamorph software by Gataca Systems. Twenty-four hours after transfection, cells were kept at 37 • C and were imaged before any activation with blue light over 17 z-stack (0.5 µm) at 561 nm (0.134 mW). The same cells were imaged at the end of the activation routine using the same 17 z-stacks while keeping the same focus. Optogenetic activations were performed every two minutes for a total duration of one hour, using the laser blue light at 491 nm (0.506 mW). We selected cells for further image quanti cation based on their visible viability, on their presence in the eld of view at the end of the experiment (some cells escaped the eld of view after one hour), and on the absence of pixels saturation (very bright, saturated clusters could appear over the time course of the experiment). All laser settings and parameters of the camera (time of exposition, gain) were kept constant for all experiments and calibration of the concentration.

A.2 Image Analysis

The initial concentration of the protein is obtained from the cell image at the initial time, t = 0. The cell is separated from the background and the intensity is computed as the average of the intensity in the cell after subtracting the background intensity, using Matlab [MATLAB, 2018].

We estimated the volume of the cells by measuring the area of the cell just above the coverslide and assuming an e ective height such that the total intensity of the 3D nal image equals the total intensity of this 2D initial image times this e ective height. This e ective height parameter varies from cell to cell and has a mean value of 1.1µm and a standard deviation of 0.4µm.

In order to quantify the size and frequency of the cluster at t = 1 hr, we smoothed the image with a gaussian lter, substracted the mean background intensity, located the local maxima of intensity in the image, and performed a watershed transform to estimate the spatial extent of each cluster [Meyer, 1994]. The size of the clusters is then determined by considering that the cluster is composed of the pixels that have at least one fth of the intensity of the maximum of such cluster. In addition, we considered a bright spot to be a cluster only if the intensity of its peak is at least 2000 arbitrary units above the background intensity-which corresponds to peaks with at least ∼ 20 monomers. Once the clusters are located and their boundaries de ned, we A.3. P F 115 add up the total intensity of each of them, separately, to obtain an estimate of the mass of each cluster, i.e., the total number of monomers in each of the clusters. This number might be slightly underestimated due to a potential self-quenching e ect of the uorescent tag upon aggregation.

A.3 Parameter Fitting

The two datasets that we want to t with Eqs. (1.28) and (1.29), i.e., cluster density and mean cluster mass, have di erent units and numerical values. In what follows, we will introduce a leastsquare minimization such that, when minimising the squares to nd the best tting parameters, both datasets are equally taken into account. To achieve this, we introduce

where f 1,2 are de ned by Eqs. (1.28) and (1.29), the 2-tuples (x i , y i ) (1,2) denote each of the datapoints i of each dataset (1 or 2, cluster concentration or cluster size), and µ 1,2 are the mean values of the datapoints of each dataset:

, M 1,2 being the number of datapoints.

A

A. A P

Supplementary Information and Figures for the nucleoid model B.1 Equilibrium free energy of an interacting system: the virial expansion

In this Appendix we brie y explain how to obtain the equilibrium free energy of a binary gas of hard spheres, as the one presented in Section 2.2.1 that yields Eq. ( 2.4). The model is composed of N A and N B particles of species A and B, which are hard spheres with radii R A , R B , respectively, con ned in a volume V. We assume that the system has cylindrical symmetry, and thus V is e ectively one dimensional, so as to apply our ndings to the model de ned by Eqs. (2.1) to (2.3).

Without taking into account reaction processes, we work out the free energy from the partition function

where x and p are the positions and momenta of the particles, Z ideal is the partition function of the system in the absence of excluded-volume interactions, i.e., an ideal gas [Huang, 1987]. In order to compute this partition function we will resort to the virial expansion [Onnes, 1902]. We introduce the Mayer function f :

where V ij is the interaction potential between the ith and jth particle. Hence, the interaction term can be expanded as follows B.3) where, assuming low density, we expect that, on average, f (x i , x j ) is small compared to one. The leading contribution in Eq. (B.3) is O(1)-the ideal gas-and the subleading one is O(f (x i , x j )), which yields the second term in the virial expansion. Here we will also compute the third virial term, O(f (x i , x j ) 2 ), because the density in our system is not low enough to allow us to safely neglect three-particle interactions. As shown in However, unlike the scaling of Section 2.3.2 Singlechromosome lamentous growth, here we consider cells which grow lamentously with a single chromosome, with constant mRNA and ribosome concentrations, independent of cell length. The predicted nucleoid length is small compared to the experimentally measured length, even for small cells.

Appendix C

Supplementary Notes on non-ideal Chemical Reaction Networks C.1 Form of the propensity functions for a regular solution theory

We now consider a model of a solution based on a lattice where each chemical species (including the solvent) occupies one lattice site, e ectively neglecting any di erences in molecular volumes. Throughout most of this manuscript we will consider a free energy with virial coe cients to correct for interactions among particles, where the non-ideality of the systems stems from.

In a lattice of size V (number of lattice sites) occupied by N di erent species, provided that N i=1 n i = V , the con gurational entropy is given by

where Ω is the number of microstates and k B is the Boltzmann constant. The internal free energy of each species is given, as before, by the standard-state chemical potential µ 0 α . Moreover, in this regular solution model, we take into account interactions among neighbouring sites, whose energy, in a mean-eld approximation with respect to space, is

where χ α,β is the energy of the interaction between a particle α and a particle β but can also be interpreted as a matrix of virial coe cients.

Taken into account the previous considerations, the free energy for a homogeneous mixture of chemical species in the regular solution model reads

With this expression of the free energy we can now derive an expression for the propensity functions, as given in eq. (3.12). The forward (or backward) rates are a function of the free energy 127