Emergence of Simple Characteristics for Heterogeneous Complex Social Agents

Abstract

Models of interacting social agents often represent agents as very simple entities with a small number of degrees of freedom, as exemplified by binary opinion models for instance. Understanding how such simple individual characteristics may emerge from potentially much more complex agents is thus a natural question. It has been proposed recently in [E. Bertin, P. Jensen, C. R. Phys. 20, 329 (2019)] that some types of interactions among agents with many internal degrees of freedom may lead to a ‘simplification’ of agents, which are then effectively described by a small number of internal degrees of freedom. Here, we generalize the model to account for agent intrinsic heterogeneity. We find two different simplification regimes, one dominated by interactions, where agents become simple and identical as in the homogeneous model, and one where agents remain strongly heterogeneous although effectively with simple characteristics.

1. Introduction

When modeling complex systems, statistical physicists often posit that the interacting entities they consider have simple individual properties, and that the possibly more complex behavior observed at a collective level, when considering many such simple entities, is simply the result of interactions between entities [1,2,3,4]. The emergence of macroscopically different properties from simpler interacting entities, often through a collective symmetry breaking mechanism, has been emphasized long ago by P. W. Anderson in his seminal paper “More is Different” [5] as one of the key mechanisms at play to account for the wealth of different objects or behaviors found in the real world. In this paper, Anderson actually presents the statistical physics approach as a generic and somewhat iterative method, to deal with the emerging complexity in a multilevel way, the outcome of a given level of analysis being the building blocks (i.e., the interacting entities) of the next level. The key role played by symmetry breaking phenomena in the emerging complexity at collective scale has been widely acknowledged in many different contexts (see, e.g., [2,6,7,8] among many others). Interestingly, a perhaps less explicit suggestion of Anderson’s paper is also to consider statistical physics models where the interacting entities are not simple objects with a handful of characteristics, but are already rather complex objects with many internal degrees of freedom. Although it has been formulated almost 50 years ago, this suggestion looks very timely by now.

In recent decades, statistical physics has gone beyond the equilibrium paradigm based on molecular entities, and has indeed started to consider assemblies of macroscopic and potentially complex objects like grains of sand [9,10], or active particles modeling for instance some types of bacteria, or self-phoretic colloids [8]. However, in these examples, macroscopic particles may in practice still be considered to be simple particles, as their many internal degrees of freedom can be subsumed into a small number of effective parameters encoding their macroscopic, non-equilibrium character, like dissipation coefficients [9,10] or self-propulsion forces [8]. To find genuine examples of assemblies of complex entities in a theoretical context, one may rather turn to the field of population dynamics and evolution in theoretical biology, where for instance large populations of individuals characterized by a complicated genome evolve under some evolutionary rules [11,12].

When considering models of social systems, it would be necessary at first sight to take into account the intrinsic complexity of human beings [1]. However, this complexity is too extreme to be captured by any type of statistical physics models, so that statistical physicists have often considered very simple models of social agents, retaining a small number of characteristics, for instance a binary [13] or continuous [14] opinion. A natural question is then to understand how more complex agents could in some situations reduce their intrinsic complexity to effectively appear as simple. This question can already be addressed within the framework of statistical physics, because there is no need to model the full complexity of human beings to address at least some aspects of this issue. A simple tentative answer has been given in [15], by illustrating on a toy model how the simplification of agents with many internal degrees of freedom may result from interactions among agents. This point of view is qualitatively consistent with the idea advocated by some sociologists that ‘the whole is less than the parts’ [16,17], in the sense that, roughly speaking, human beings may leave aside part of their complexity to build a group. It has been argued in [15] that this point of view can be reconciled with the seemingly antagonist viewpoint of statistical physicists according to which ‘the whole is more than the sum of the part’, due to collective phenomena and symmetry breakings.

Yet, the model introduced in [15] considered an assembly of identical agents, while agent heterogeneity may be expected to be an important characteristic of human beings. In this note, we extend the model of [15] by including heterogeneity between agents. We show, using methods inspired by the physics of glasses, that two different types of agents simplification can occur in this model, one driven by interactions as in [15], and the other one driven by heterogeneity.

2. Model

We introduce a model of complex agents generalizing the one introduced in [15] by now considering heterogeneous agents. The model is composed of N interacting agents with an internal state described by a configuration ${\mathcal{C}}_i\in \{1,\dots ,H\}$, with $i=1,\dots ,N$ and where the number H of configurations is large. We write for later convenience H in the form $H=n_0^M$, where $n_0$ is a fixed integer number, and M is assumed to be large. This is typically the case if the configuration ${\mathcal{C}}_i$ is composed of M degrees of freedom, each of which taking $n_0$ possible values. Each agent is endowed with a characteristic that can be either present or absent depending on the configuration ${\mathcal{C}}_i$. Intuitively, this characteristic could be a preference for a specific kind of music or tempo for the members of a vocal ensemble [17] for instance, or more generally be related to a binary opinion [13]. This feature is conveniently encoded by a variable $S_i\left({\mathcal{C}}_i\right)\in \{0,1\}$: the characteristic is present when $S_i\left({\mathcal{C}}_i\right)=1$, and absent when $S_i\left({\mathcal{C}}_i\right)=0$. The general idea is that the characteristic would be present only in a small number of internal states, so that it would typically remain unobserved except if there is a strong probability bias towards the few configurations for which $S_i\left({\mathcal{C}}_i\right)=1$. To implement this idea in practice in the model, we assume that the characteristic is present in a single configuration, which we label ${\mathcal{C}}_i=1$.

We now need to define the dynamics of agents. Following standard practice in the modeling of social agents [18,19], we assume that the dynamics is driven by an individual utility function $u_i=u_i\left({\mathcal{C}}_i\right|{\mathcal{C}}_{j\ne i})$ that accounts both for individual preferences and for interactions with other agents. An agent i stochastically changes configurations according to the following rule. Given the current configuration ${\mathcal{C}}_i$, the new configuration ${\mathcal{C}}_i^{\prime }$ is randomly chosen with a probability rate given by the logit rule, a standard rule in the field of social agents modeling [18,19],

$$W\left({\mathcal{C}}^{\prime }\right|\mathcal{C})=\frac{1}{1+e^{-\mathsf{\Delta }{u}_i/T}}$$

(1)

with $\mathsf{\Delta }{u}_i=u_i\left({\mathcal{C}}_i^{\prime }\right)-u_i\left({\mathcal{C}}_i\right)$ the variation of utility generated by the change of configuration (note that configurations ${\mathcal{C}}_j$ of the other agents $j\ne i$ are kept fixed). The parameter T plays a role similar to temperature in statistical physics, and characterizes the degree of stochasticity in the decision rule.

Our goal is to model heterogeneous agents that interact through their characteristic $S_i$ (their internal state is otherwise invisible to other agents). Such an interaction between characteristics was already introduced in Ref. [15] for the homogeneous case, by considering a mean-field (or fully connected) coupling. We keep here the same mean-field expression of the interaction contribution to the utility. To model heterogeneity in a minimal way, we include an additive contribution $U_i\left({\mathcal{C}}_i\right)$ to the utility of agent i, which corresponds to the agent’s intrinsic (or idiosyncratic) utility associated with configuration ${\mathcal{C}}_i$. Under these assumptions, we end up with the following form of the utility function of agent i,

$$u_i\left({\mathcal{C}}_i\right|{\mathcal{C}}_{j\ne i})=U_i\left({\mathcal{C}}_i\right)+\frac{K}{N}\sum _{j(\ne i)}{S}_i\left({\mathcal{C}}_i\right)S_j\left({\mathcal{C}}_j\right),$$

(2)

where K is the coupling constant characterizing the interaction with the other agents. Note the $1/N$ scaling of the interaction term, typical of fully connected models where all particles or agents interact with each other in a similar way. Another, and possibly interesting way to introduce heterogeneity would be to have heterogeneous couplings among agents, by replacing the single coupling constant J by different coupling constants $J_{ij}$ that would depend on the pair of agents $(i,j)$. However, this would make the analysis of the model much more involved. In the present work, we restrict ourselves to the simplest case when heterogeneity is introduced only in the idiosyncratic utility.

To fully define the model, we need to specify the heterogeneous intrinsic utilities $U_i\left({\mathcal{C}}_i\right)$. Inspired by spin-glass physics, and more precisely by the Random Energy Model for spin-glasses [20,21], we model the heterogeneous intrinsic utilities $U_i\left({\mathcal{C}}_i\right)$ as a quenched randomness. To be specific, for all $i=1,\dots ,N$ and ${\mathcal{C}}_i=1,\dots ,H$, the intrinsic utility $U_i\left({\mathcal{C}}_i\right)$ is randomly drawn from a Gaussian distribution $\rho \left(U\right)$,

$$\rho \left(U\right)=\frac{1}{\sqrt{\pi M{J}^2}}{e}^{-U^2/M{J}^2},$$

(3)

where we recall that M is defined by $H=n_0^M$. The utilities $U_i\left({\mathcal{C}}_i\right)$ do not change in time. The Gaussian form of the intrinsic utility distribution, as well as the scaling of the variance with M, is guided by similar assumptions made on the energy distribution in the original version of the Random Energy Model [20,21]. However, the results obtained are not expected to depend on this Gaussian assumption at a qualitative level, because the analysis conducted here is similar to the one performed in the Random Energy Model, whose phenomenology is the same for a broad class of energy distributions [22,23].

Now that the dynamics of the model has been defined, we would like to determine the steady-state probability distribution $P({\mathcal{C}}_1,\dots ,{\mathcal{C}}_N)$ of the model. Although the transition rate (1) looks similar to the Glauber or Kawasaki rates used in statistical physics for spin models or lattice gas models, it is actually different because the utility $u_i\left({\mathcal{C}}_i\right|{\mathcal{C}}_{j\ne i})$ is a function attached to each agent i, and not a unique global function like the energy in physics models (see, e.g., [24] for a more detailed discussion of this issue). Yet, it is possible in some cases to find a global function $E({\mathcal{C}}_1,\dots ,{\mathcal{C}}_N)$ that thus does not depend on the agent i, such that the variation $\mathsf{\Delta }E$ during a transition from one configuration to another is equal to the variation $\mathsf{\Delta }{u}_i$ of the utility of agent i. The global observable E plays a role similar to the energy in physics (up to a change of sign). The main difference is that in physics models, the energy (or Hamiltonian) is given by the definition of the model, whereas in socioeconomic models the function E has to be determined (if it exists) from the dynamical rules defined in terms of the utility. Here, due to the simple enough form of the utility function, such a global function E exists and can be determined explicitly as

$$E({\mathcal{C}}_1,\dots ,{\mathcal{C}}_N)=\sum _{i=1}^N{U}_i\left({\mathcal{C}}_i\right)+\frac{K}{2N}\sum _{i,j(i\ne j)}{S}_i\left({\mathcal{C}}_i\right)S_j\left({\mathcal{C}}_j\right).$$

(4)

Please note that the quantity E is thus different from the total utility ${\sum }_i{u}_i$, due to the factor $\frac{1}{2}$ in the interaction term.

Given the property $\mathsf{\Delta }{u}_i=\mathsf{\Delta }E$, it can now be easily checked that the dynamics defined by the transition rate Equation (1) obeys the detailed balance property in terms of the equilibrium distribution

$$P({\mathcal{C}}_1,\dots ,{\mathcal{C}}_N)=\frac{1}{Z}{e}^{\beta E({\mathcal{C}}_1,\dots ,{\mathcal{C}}_N)},$$

(5)

with $\beta \equiv 1/T$, and where Z is a normalization constant. In this form, it becomes clear that the present model shares similarities with the Random Energy Model [20,21,22] as mentioned above, as well as with the Ising model [7], the Potts model [25] or the Blume-Emery-Griffiths spin-1 model [26]. However, it also exhibits important differences with each of these models.

Now that we determined the equilibrium distribution of the model, our goal is to investigate its phase diagram to assess the effect of the competition between agents heterogeneity (due to their quenched intrinsic utility) and collective effects that could arise from interactions. As recalled in the introduction, the homogeneous version of the model, studied in [15], exhibits a transition driven by interactions between a phase where agents essentially visit all their internal states and thus have no strongly preferred configurations, and an ordered state where all agents ‘standardize’ in the same configuration ${\mathcal{C}}_1$ such that $S_i\left({\mathcal{C}}_1\right)=1$ (i.e., the characteristic of the agent becomes visible). In the present generalization of the model, we wish to explore whether the standardized state survives the heterogeneity of agent intrinsic utility. Indeed, the internal state ${\mathcal{C}}_1$ such that $S_i=1$ may have a lower intrinsic utility $U_i\left({\mathcal{C}}_1\right)$ than other more favored configurations of agent i, and this may prevent a common standardization of all agents in the same state.

To investigate this issue, we introduce the order parameter q defined as

$$q=\frac{1}{N}\sum _{i=1}^N{S}_i.$$

(6)

If the number N of agents is large, the pseudo-energy E can be expressed in terms of the order parameter q as

$$E=\sum _{i=1}^N{U}_i+\frac{1}{2}NK{q}^2.$$

(7)

We now wish to determine the distribution $P\left(q\right)$ of the order parameter q. This distribution is obtained by summing the joint distribution $P({\mathcal{C}}_1,\dots ,{\mathcal{C}}_N)$ over all configurations $({\mathcal{C}}_1,\dots ,{\mathcal{C}}_N)$ with a given value of q, namely

$$P\left(q\right)=\sum _{{\mathcal{C}}_1,\dots ,{\mathcal{C}}_N}P({\mathcal{C}}_1,\dots ,{\mathcal{C}}_N)\delta \left(q({\mathcal{C}}_1,\dots ,{\mathcal{C}}_N)-q\right)$$

(8)

with $\delta$ a generalized Kronecker delta such that $\delta \left(0\right)=1$ and $\delta \left(x\right)=0$ otherwise. For a given value of q, there are $qN$ agents in configuration ${\mathcal{C}}_i=1$ and $(1-q)N$ agents in any of the other $n_0^M-1$ configurations. To evaluate $P\left(q\right)$ according to its definition given in Equation (8), we thus need to sum over all possible ways to choose the $qN$ agents in configuration ${\mathcal{C}}_i=1$ and the $(1-q)N$ agents in other configurations. Technically, we need to sum over all possible unordered subsets with $qN$ elements of the integer set $[1,N]$. To write reasonably compact equations, we denote by ${\mathcal{S}}_{q,N}$ a generic unordered subset with $qN$ elements of the set $[1,N]$, and we sum over all possible such subsets ${\mathcal{S}}_{q,N}$. With these notations, $P\left(q\right)$ can be written as

$$P\left(q\right)=\frac{1}{Z}\sum _{{\mathcal{S}}_{q,N}}\left[\prod _{i\in {\mathcal{S}}_{q,N}}{e}^{\beta U_i\left(1\right)}\right]\left[\prod _{i\notin {\mathcal{S}}_{q,N}}\sum _{{\mathcal{C}}_i=2}^{n_0^M}{e}^{\beta U_i\left({\mathcal{C}}_i\right)}\right]e^{\frac{1}{2}\beta NK{q}^2},$$

(9)

where in the second product, the index i is implicitly restricted to the interval $[1,N]$. The distribution $P\left(q\right)$ defined in Equation (8) actually depends on the specific realization of the random utilities $U_i\left({\mathcal{C}}_i\right)$, so that $P\left(q\right)$ should in principle be eventually averaged over these random utilities. However, to make calculations easier, we do not compute explicitly the average over $U_i\left({\mathcal{C}}_i\right)$, but rather use heuristic arguments to evaluate the typical values of the random quantities appearing in Equation (9). Such an estimate is expected to be sufficient to determine the leading exponential behavior $P\left(q\right)\sim e^{-Nf\left(q\right)}$ of the distribution $P\left(q\right)$. The first product between brackets in Equation (9) can be rewritten as the exponential of ${\sum }_{i\in {\mathcal{S}}_{q,N}}\beta U_i\left(1\right)$. The latter sum is (up to the factor $\beta$) a sum of independent and identically distributed Gaussian random variables drawn from the distribution Equation (3). The sum is thus also Gaussian distributed, with zero mean and variance $\frac{1}{2}qNm{\left(\beta J\right)}^2$. It follows that

$$\prod _{i\in {\mathcal{S}}_{q,N}}{e}^{\beta U_i\left(1\right)}\sim e^{O\left(\sqrt{N}\right)},$$

(10)

and this product can thus be neglected (i.e., replaced by 1) when looking for the behavior of $P\left(q\right)$ at exponential order in N.

The key point to determine $P\left(q\right)$ explicitly is now to evaluate the typical value ${\mathcal{Z}}_{\mathrm{typ}}$ of the sum

$${\mathcal{Z}}_i=\sum _{{\mathcal{C}}_i=2}^{n_0^M}{e}^{\beta U_i\left({\mathcal{C}}_i\right)},$$

(11)

where we recall that each $U_i\left({\mathcal{C}}_i\right)$ is an independent random variable drawn from the distribution $\rho \left(U\right)$ given in Equation (3). Once this estimate is known, the distribution $P\left(q\right)$ can be approximated as

$$P\left(q\right)\sim \frac{1}{Z}\left(\genfrac{}{}{0pt}{}{N}{qN}\right){\mathcal{Z}}_{\mathrm{typ}}^{(1-q)N}{e}^{\frac{1}{2}\beta NK{q}^2}.$$

(12)

The quantity ${\mathcal{Z}}_i$ has the same form as the partition function of the Random Energy Model (REM) [20,21], and we can thus borrow some methods from the REM to evaluate it. A standard approach to study the REM is to evaluate the density of states of a typical realization of the disorder. It is convenient at this stage to define the rescaled variable $u=U/M$. Denoting as $n\left(u\right)$ the density of state of a given realization, the average $\langle n\left(u\right)\rangle$ over the disorder is given by

$$\langle n\left(u\right)\rangle =n_0^M\tilde{\rho }\left(u\right)\sim e^{M(ln{n}_0-u^2/J^2)}$$

(13)

to exponential order in M, with the definition $\tilde{\rho }\left(u\right)=M\rho \left(Mu\right)$. Hence if $ln{n}_0-u^2/J^2>0$, corresponding to $\left|u\right|<u_0=J\sqrt{ln{n}_0}$, the average density of state $\langle n\left(u\right)\rangle$ is exponentially large with M, and for a typical sample $n\left(u\right)\approx \langle n\left(u\right)\rangle$. By contrast, when $\left|u\right|>u_0$, the average density of state is exponentially small with M, meaning that in most realizations there are actually no states with $\left|u\right|>u_0$ (the exponentially small value of the average density of states comes from very rare and atypical realizations with a few states in this range). Therefore, to describe a typical realization, one can in practice consider that the density of state is equal to zero for $\left|u\right|>u_0$, and equal to $\langle n\left(u\right)\rangle$ for $\left|u\right|<u_0$. We can thus evaluate ${\mathcal{Z}}_i$ as

$${\mathcal{Z}}_i=\sum _{{\mathcal{C}}_i}{e}^{\beta U_i\left({\mathcal{C}}_i\right)}\approx \sqrt{\frac{M}{\pi J}}{\int }_{-u_0}^{u_0}du{e}^{M(ln{n}_0-u^2/J^2+\beta u)}\equiv {\mathcal{Z}}_{\mathrm{typ}}.$$

(14)

The integral in Equation (14) can be evaluated by a saddle-point calculation in the large M limit. Defining $g\left(u\right)=ln{n}_0-u^2/J^2+\beta u$, ${\mathcal{Z}}_{\mathrm{typ}}$ is given for large M by ${\mathcal{Z}}_{\mathrm{typ}}\sim e^{M{g}_{max}}$, where $g_{max}$ is the maximum value of $g\left(u\right)$ over the interval $[-u_0,u_0]$. Let us first look for the maximum of $g\left(u\right)$ over the whole real axis. Defining $u^{*}$ such that $g^{\prime }\left(u^{*}\right)=0$, we find $u^{*}=\frac{1}{2}\beta J^2$. Recalling that $\beta =1/T$, and defining

$$T_{\mathrm{g}}=\frac{J}{2\sqrt{ln{n}_0}},$$

(15)

we find that $-u_0<u^{*}<u_0$ for $T>T_{\mathrm{g}}$, whereas $u^{*}\ge u_0$ for $T\le T_{\mathrm{g}}$. Hence $g_{max}=g\left(u^{*}\right)$ for $T>T_{\mathrm{g}}$, and $g_{max}=g\left(u_0\right)$ for $T\le T_{\mathrm{g}}$, taking into account the fact that $g\left(u\right)$ is monotonously increasing for $u<u^{*}$. We thus obtain

$${\mathcal{Z}}_{\mathrm{typ}}\approx \left\{\begin{array}{cc}{e}^{M(ln{n}_0+\frac{1}{4}{\beta }^2{J}^2)}\hfill & \mathrm{if}T>T_{\mathrm{g}},\hfill \\ e^{M\beta J\sqrt{ln{n}_0}}\hfill & \mathrm{if}T\le T_{\mathrm{g}}.\hfill \end{array}\right.$$

(16)

Using Equations (12) and (16), we obtain after expanding the factorials thanks to the Stirling formula that $P\left(q\right)$ takes a large deviation form $P\left(q\right)\sim e^{Nf\left(q\right)}$. In physical terms, $f\left(q\right)$ may be thought of as (the opposite of) a free energy. The explicit expression of $f\left(q\right)$ depends on the temperature range. For $T>T_{\mathrm{g}}$, $f\left(q\right)$ reads

$$f\left(q\right)=-qlnq-(1-q)ln(1-q)+(1-q)M\left[ln{n}_0+\frac{J^2}{4T^2}\right]+\frac{1}{2T}K{q}^2,$$

(17)

whereas for $T\le T_{\mathrm{g}}$,

$$f\left(q\right)=-qlnq-(1-q)ln(1-q)+(1-q)M\frac{J}{T}\sqrt{ln{n}_0}+\frac{1}{2T}K{q}^2.$$

(18)

Please note that we did not take into account in $f\left(q\right)$ the contribution coming from the normalization factor Z, as this would simply add a constant to $f\left(q\right)$. Inspired by Ref. [15], where interesting results were obtained for a coupling constant $K\sim M$, we assume in what follows that

$$K=kM,$$

(19)

and take the reduced constant k as the relevant control parameter in the model (on top of temperature T). For large M, the expression of $f\left(q\right)$ then simplifies to

$$f\left(q\right)=\left\{\begin{array}{cc}\frac{M}{T}\left[(1-q)\left(Tln{n}_0+\frac{J^2}{4T}\right)+\frac{1}{2}k{q}^2\right]\hfill & \mathrm{if}T>T_{\mathrm{g}},\hfill \\ \frac{M}{T}\left[(1-q)J\sqrt{ln{n}_0}+\frac{1}{2}k{q}^2\right]\hfill & \mathrm{if}T\le T_{\mathrm{g}}.\hfill \end{array}\right.$$

(20)

We first observe that in this large-M approximation, $f\left(q\right)$ is a convex function of q for all values of temperature T, so that the maximum of $f\left(q\right)$ over the interval $0\le q\le 1$ is either $f\left(0\right)$ or $f\left(1\right)$. The most probable state is found to be $q=0$ for $k<k_{\mathrm{c}}\left(T\right)$ and $q=1$ for $k>k_{\mathrm{c}}\left(T\right)$, where the critical line $k_{\mathrm{c}}\left(T\right)$ is defined as (with $k_0=2J\sqrt{ln{n}_0}$)

$$k_{\mathrm{c}}\left(T\right)=\left\{\begin{array}{cc}\frac{k_0}{2}\left(\frac{T}{T_{\mathrm{g}}}+\frac{T_{\mathrm{g}}}{T}\right)\hfill & \mathrm{if}T>T_{\mathrm{g}},\hfill \\ k_0\hfill & \mathrm{if}T\le T_{\mathrm{g}}.\hfill \end{array}\right.$$

(21)

The curve $k_{\mathrm{c}}\left(T\right)$ thus separates the $(k,T)$ phase diagram into two regions, a region with $q=0$ at low coupling and a region with $q=1$ at high coupling. The corresponding phase diagram is plotted in Figure 1. Please note that for $J=0$ (i.e., in the absence of disorder), the glassy region in the phase diagram disappears, and one recovers the phase transition at $T_{\mathrm{c}}=k/(2ln{n}_0)$ between a high-temperature phase with $q=0$ and a low-temperature phase with $q=1$ found in the homogeneous model [15].

Moreover, we have also seen that a change of behavior occurs at $T=T_{\mathrm{g}}$. For $T>T_{\mathrm{g}}$, the agents dynamically visit a large number of configurations, while for $T<T_{\mathrm{g}}$ their dynamics becomes essentially frozen, and only a few configurations have a significant probability to be visited. In other words, for $T<T_{\mathrm{g}}$ agents become ‘stuck’ in a small number of configurations with the highest utility. In the context of the Random Energy Model for glasses, the temperature $T_{\mathrm{g}}$ corresponds to the glass transition.

Hence there are actually three different regions in the phase diagram shown in Figure 1. For $T>T_{\mathrm{g}}$ and $k<k_{\mathrm{c}}\left(T\right)$, agents have no preferred configurations and visit many different configurations over time. For $T<T_{\mathrm{g}}$ and $k<k_{\mathrm{c}}\left(T_{\mathrm{g}}\right)$, each agent spends a lot of time in a small set of preferred configurations. In other words, agents look simple, but they remain different one from the other. This regime is dominated by agent heterogeneity, and there is on average no macroscopic overlap between agent configurations ($q=0$). In the last region $k>k_{\mathrm{c}}\left(T\right)$, the coupling between agents dominates over agents heterogeneity, and all agents are essentially in the same configuration, leading to a strong overlap ($q=1$) and to the emergence of a common characteristic, a phenomenon that has been called standardization in [15].

3. Conclusions

We have extended here the model proposed in [15], where agents have many internal configurations, but can select a specific configuration thanks to interactions, leading to simplified, or standardized agents. Although Ref. [15] focused on homogeneous agents, we have in the present work extended the model to account for agent heterogeneity, introduced through random idiosyncratic utilities associated with each configurations of each agent. Heterogeneity is introduced in a minimal way, directly inspired from the Random Energy Model for glasses [20]. Including heterogeneity in the model leads to the onset of a new phase at low temperature and small coupling, when heterogeneity dominates over interactions among agents. In physical terms, this new phase shares similarity with a glass phase. The main difference here is that we consider not a single glassy system as in physics, but a large assembly of interacting agents, each of which experiencing an internal glass transition. In this glassy phase, agents are stuck in a small number of possible internal configurations, those with maximal idiosyncratic utility. Despite the individual simplification of agents, the assembly remains heterogeneous, because all agents are in different configurations and do not share common characteristics. Increasing the coupling strength, interactions eventually dominate over heterogeneity, effectively leading to simple agents all occupying the same configuration and sharing the same characteristic. The latter phase is similar to the low-temperature phase found in the homogeneous model of [15], where agents have been called standardized. Finally, at high temperature, agents keep their internal complexity and do not significantly feel interactions or heterogeneity, again similarly to the homogeneous model [15].

A further step in the study of this model would be to study possible collective effects that could arise once agents are simplified. This has been done in the homogeneous model by assuming that standardized agents could be in one of two distinct states, corresponding for instance to two different opinions (or to a spin in physical terms). Interactions between these binary degrees of freedom may then lead, at low temperature, to a collectively ordered state. It would be interesting to see in more details how agent heterogeneity could possibly modify this simple picture. In addition, the role of a lower connectivity (in the present work, agents interact with all other agents) would clearly deserve to be investigated, as connectivity is known in many situations to modify critical properties [7].