# Exact solutions of the simplified March model for organizational learning

Hang-Hyun Jo<EMAIL_ADDRESS>Department of Physics, The Catholic
University of Korea, Bucheon 14662, Republic of Korea

###### Abstract

James G. March’s celebrated agent-based simulation model for organizational
learning [March, Organization Science 2, 71 (1991)] has been extensively
studied for the last decades. Yet the model was not fully understood due to
the lack of analytical solutions of the model. We simplify the March model to
take an analytical approach using master equations. We then derive exact
solutions for some simplest yet nontrivial cases, and perform numerical
estimation of master equations for more complicated cases. Both analytical and
numerical results are in good agreement with agent-based simulations. These
results are also compared to those of the original March model. Our approach
enables us to rigorously understand the results of the simplified model as
well as the original model to a large extent.

††preprint: APS/123-QED

##  I Introduction

James G. March introduced an agent-based simulation model for organizational
learning in his seminal paper in 1991 [1]. Since then, the original March
model, its variants, and other similar models have been extensively studied
for the last decades [2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15]. March’s
model considers an external reality, an organizational code (code hereafter),
and individual members of the organization. The code represents a set of
norms, rules, etc. that is updated using the knowledge of individuals about
the reality, while individuals learn from the code about the reality too. By
doing so, the organizational knowledge about the reality is collected from
individuals and it is disseminated to them at the same time. He studied the
effects of learning rates of individuals and of the code on their achieved
knowledge about the reality. In order to consider more realistic situations,
he took personnel turnover and environmental turbulence into account in his
model to find that there may exist an optimal turnover rate that maximizes the
achieved knowledge, depending on the learning rates.

We remark that the March model can be considered in the framework of opinion
dynamics in networks [16, 17, 18, 19, 20, 21, 22, 23, 24]. That is, the code
plays a hub node in a hub-and-spoke network, while individuals are dangling
nodes [25]. Nodes update their opinions or beliefs according to their
neighbors’ opinions or beliefs, while the external reality plays an external
field or source affecting all nodes. It implies that various analytical
approaches developed for the opinion dynamics can be applied to the March
model.

The March model and its variants have provided with insights into management
and business administration [8] but mostly by means of computer simulations of
models [14]. In general, for rigorous understanding of models, derivation of
their exact, analytical solutions is of utmost importance. In our work, we
simplify March’s original model to explicitly write down master equations
describing the dynamics of the model. Then we derive exact solutions of the
simplified model for some simplest yet nontrivial cases. Numerical estimation
of master equations is performed for more complicated cases. Both analytical
and numerical results are shown to be in good agreement with agent-based
simulation results. Our approach enables us to rigorously understand the
results of not only the simplified model but also the original model to a
large extent.

The paper is organized as follows. In Sec. II, we describe the original March
model and our simplified version. In Sec. III, analytical, numerical, and
simulation results of the simplified model are presented and compared to the
results of the original March model. Finally, we conclude our work in Sec. IV.

##  II Models

###  II.1 Original March model

As mentioned, the original model by March considers the external reality, the
code, and individual members of the organization [1]. We remark that in this
Subsection we will use mathematical symbols originally used in the March’s
paper, and they are not supposed to be confused with symbols in the next
Subsection and throughout the paper. The model is based on the following
assumptions.

(i) The reality is characterized in terms of an $m$-dimensional vector, each
element of which may have the value of $1$ or $-1$ with equal probabilities.

(ii) The code and $n$ individuals in the organization have beliefs about the
reality. The belief is also represented by an $m$-dimensional vector, each
element of which may have the value of $1$, $0$, or $-1$ with equal
probabilities. These beliefs may change over time.

(iii) At each time step, each individual may change elements of its belief
that are different from those of the code unless the code’s element is $0$.
Each of such elements of the individual changes to that of the code with a
probability $p_{1}$, independently of other elements.

(iv) At the same time, the code updates its belief based on beliefs of some
individuals. For this, individuals whose beliefs are closer to the reality
than the code is to the reality are identified, which are called a superior
group. Then each element of the code’s belief changes to the dominant element
within the superior group with a probability $p_{2}$, independently of other
elements.

So far, the reality has been assumed to be fixed and the individuals are not
replaced by new ones. Thus, it is called a closed system. March first
considered a homogeneous population in the closed system, in which all
individuals are assigned the identical learning probability. Then he
considered the heterogeneous population in the closed system, such that some
individuals have higher learning probability than the others. Finally, a
homogeneous population in the open system is also considered; in the open
system individuals may be replaced by new ones (turnover) and/or the reality
changes over time (turbulence). The turnover probability is denoted by $p_{3}$
and the turbulence probability is by $p_{4}$. That is, with a probability
$p_{3}$ each individual is replaced by a new one having a random belief vector
at each time step. Also, each element of the reality shifts to the other
value, i.e., from $1$ to $-1$ or from $-1$ to $1$, with a probability $p_{4}$.

###  II.2 Simplified March model

Let us simplify the March model. As in the original model we consider an
external reality, a code, and $N$ agents. At a time step $t$, the external
reality, denoted by a variable $r(t)$, can have a value of $0$ or $1$. Beliefs
of the code and agents about reality are respectively represented by
variables, namely, $c(t)\in\\{0,1\\}$ for the code and
$\sigma_{i}(t)\in\\{0,1\\}$ for the $i$th agent with $i=1,\ldots,N$. For a
given initial condition of $r(0)$, $c(0)$, and $\\{\sigma_{i}(0)\\}$, each
time step consists of four stages:

(i) Every agent of $i\in\\{1,\ldots,N\\}$ independently updates their belief
by learning from the code with a _socialization_ probability
$p^{(i)}\in(0,1]$:

$\displaystyle\sigma_{i}(t+1)=c(t).$ (1)

(ii) Each agent is replaced by a new agent with a _turnover_ probability
$u\in[0,1]$, and the new agent is assigned a belief randomly drawn from
$\\{0,1\\}$.

(iii) The code learns from agents who are superior to the code. Here superior
agents indicate those whose beliefs are closer to the reality than the code’s
belief is. For example, if $r(t)=1$, the code learns from superior agents only
when $c(t)=0$ and there is at least one agent with $\sigma_{i}(t)=1$. Denoting
a superior agent to the code by $j$, the code updates its belief with a
_codification_ probability $q\in(0,1]$:

$\displaystyle c(t+1)=\sigma_{j}(t)\ \textrm{for}\
j\in\\{i|\delta_{\sigma_{i}(t),r(t)}>\delta_{c(t),r(t)}\\},$ (2)

where $\delta_{\cdot,\cdot}$ is a Kronecker delta.

(iv) With a _turbulence_ probability $v\in[0,1]$, the reality is assigned a
new value randomly drawn from $\\{0,1\\}$, which closes the time step.

Since the reality, the code, and agents update their value or beliefs
synchronously, the order of four stages does not affect the result, except for
the socialization and turnover of agents. Note that parameters
$\\{p^{(i)}\\}$, $q$, $u$, and $v$ in our simplified model correspond to
$p_{1}$, $p_{2}$, $p_{3}$, and $p_{4}$ in the original March model [1],
respectively.

Our simplified model can be called a closed system if $u=v=0$, otherwise it is
an open system. In the open system, the reality can vary over time
(turbulence), and/or agents can be replaced by new agents (turnover). In
contrast, in the closed system, the reality is assumed to have the fixed value
of $1$ for the entire period of time, i.e., $r(t)=1$ for all $t$, without loss
of generality, and there is no turnover of agents.

##  III Results

Figure 1: (a) Schematic diagram of the homogeneous learning model in a closed
system with $r=1$, showing how the code and agents learn from each other with
learning probabilities $p$ and $q$ [Eqs. (1) and (2)]. (b) Transition
structure between states of the system. Each state is denoted by $(c,n)$,
where $c$ is the belief of the code and $n$ is the number of agents whose
belief is $1$. Self-loops are not shown for better visualization. Both $(0,0)$
and $(1,N)$ are absorbing states. (c) Analytic solutions of $\rho(\infty)$ in
Eq. (14) (lines) with simulation results (symbols) for the case with $N=40$
and an initial condition that $P_{cn}(0)=\delta_{c,0}\delta_{n,1}$. (d)
Numerical estimation of $\rho(T)$ using the master equation in Eq. (7) (lines)
for the case with $N=40$ and an initial condition that
$P_{cn}(0)=\frac{1}{2}\delta_{n,N/2}$ for each $c\in\\{0,1\\}$. Here $T$ is
the first time step satisfying $|\rho(T)-\rho(T-1)|<10^{-4}$. The
corresponding simulation results are shown with symbols. In panels (c, d),
each symbol was averaged over $10^{4}$ different runs. Standard errors are
omitted as they are smaller than symbols.

###  III.1 Homogeneous learning in a closed system

We consider the homogeneous learning model in a closed system. That is,
$r(t)=1$ for all $t$, there is no turnover of agents, and every agent has the
same learning probability $p$, i.e.,

$\displaystyle p^{(i)}=p\ \textrm{for}\ i=1,\ldots,N.$ (3)

See also Fig. 1(a). The state of the system at each time step $t$ can be
summarized in terms of the code’s belief $c(t)$ and the number of agents whose
belief matches the reality, which we denote by $n(t)\in\\{0,\ldots,N\\}$.
Precisely, $n(t)$ is defined as

$\displaystyle n(t)\equiv\sum_{i=1}^{N}\delta_{\sigma_{i}(t),r(t)}.$ (4)

In our case with $r(t)=1$, one simply has $n(t)=\sum_{i}\sigma_{i}(t)$. Then
the expected density of agents with the belief matching the reality is given
by

$\displaystyle\rho(t)\equiv\left\langle\frac{n(t)}{N}\right\rangle,$ (5)

which can be interpreted as the expected belief of a randomly chosen agent, or
average individual knowledge [3].

Depending on the initial belief of the code, two scenarios are possible.
Firstly, if $c(0)=1$, the code does not change its belief because it already
coincides with the reality, and agents’ beliefs will eventually converge to
the value of $1$ by Eq. (1). It implies an absorbing state that the code and
all agents share the same value as the reality, which is denoted by
$(c,n)=(1,N)$. Secondly, if $c(0)=0$, $n(t)$ will decrease until the code’s
belief changes to $1$ by Eq. (2) as long as there is at least one agent with
belief of $1$. Once the code’s belief becomes $1$, $n(t)$ will increase to
reach the absorbing state $(c,n)=(1,N)$. However, this is not always the case;
$n(t)$ may reach $0$ before $c(t)$ changes to $1$, implying that both the code
and agents have the belief of $0$ without further dynamics. This indicates
another absorbing state $(c,n)=(0,0)$. Figure 1(b) shows the transition
structure between states with two absorbing states emphasized in red.

For the analysis, let us denote by $P_{cn}(t)$ the probability that at time
step $t$ the code’s belief is $c$ and there are exactly $n$ agents with belief
of $1$. These probabilities satisfy the normalization condition as

$\displaystyle\sum_{c=0}^{1}\sum_{n=0}^{N}P_{cn}(t)=1.$ (6)

They evolve according to the following master equation in discrete time:

$\displaystyle
P_{cn}(t+1)=\sum_{c^{\prime}n^{\prime}}W_{c^{\prime}n^{\prime}\to
cn}P_{c^{\prime}n^{\prime}}(t),$ (7)

where the transition probabilities read [see Fig. 1(b)]

$\displaystyle\begin{split}&W_{00\to cn}=\delta_{0,c}\delta_{0,n},\\\
&W_{0n^{\prime}(\neq 0)\to 0n}=\begin{cases}{n^{\prime}\choose
n}p^{n^{\prime}-n}\bar{p}^{n}\bar{q}&\textrm{if}\ n^{\prime}\geq n,\\\
0&\textrm{if}\ n^{\prime}<n,\end{cases}\\\ &W_{0n^{\prime}(\neq 0)\to
1n}=\begin{cases}{n^{\prime}\choose
n}p^{n^{\prime}-n}\bar{p}^{n}q&\textrm{if}\ n^{\prime}\geq n,\\\
0&\textrm{if}\ n^{\prime}<n,\end{cases}\\\ &W_{1n^{\prime}\to 0n}=0\ \forall
n^{\prime},n,\\\ &W_{1n^{\prime}\to 1n}=\begin{cases}{N-n^{\prime}\choose
N-n}p^{n-n^{\prime}}\bar{p}^{N-n}&\textrm{if}\ n^{\prime}\leq n,\\\
0&\textrm{if}\ n^{\prime}>n.\end{cases}\end{split}$ (8)

Here we have used $\bar{p}\equiv 1-p$ and $\bar{q}\equiv 1-q$. Calculating Eq.
(7) recursively with any initial condition $\\{P_{cn}(0)\\}$, one can in
principle obtain $P_{cn}(t)$ for any $c$, $n$, and $t$, hence $\rho(t)$ in Eq.
(5), i.e.,

$\displaystyle\rho(t)=\frac{1}{N}\sum_{c=0}^{1}\sum_{n=0}^{N}nP_{cn}(t).$ (9)

We focus on steady states of the model. It is obvious that all initial
probabilities eventually end up with two absorbing states, i.e., $(c,n)=(0,0)$
and $(1,N)$, implying $P_{00}(\infty)+P_{1N}(\infty)=1$. Thus, the average
individual knowledge in Eq. (9) reads

$\displaystyle\rho(\infty)=P_{1N}(\infty)=1-P_{00}(\infty).$ (10)

As the simplest yet nontrivial case, let us consider an initial condition that
$P_{cn}(0)=\delta_{c,0}\delta_{n,1}$, namely, $P_{01}(0)=1$ and $P_{cn}(0)=0$
for all other states $(c,n)\neq(0,1)$. Using

$\displaystyle W_{01\to 01}=\bar{p}\bar{q}\equiv\alpha\ \textrm{and}\ W_{01\to
00}=p\bar{q}\equiv\beta,$ (11)

the master equations for $P_{01}$ and $P_{00}$ [Eq. (7)] are written as
follows:

$\displaystyle\begin{split}&P_{01}(t+1)=\alpha P_{01}(t),\\\
&P_{00}(t+1)=P_{00}(t)+\beta P_{01}(t).\end{split}$ (12)

Master equations for all other states than $P_{01}(t)$ and $P_{00}(t)$ are
irrelevant to calculate $P_{00}(\infty)$ in Eq. (10). Since
$P_{01}(t)=\alpha^{t}$, one obtains

$\displaystyle P_{00}(t)=\beta(1+\alpha+\ldots+\alpha^{t-1}),$ (13)

leading to

$\displaystyle\rho(\infty)=1-\frac{\beta}{1-\alpha}=\frac{q}{p+q-pq}.$ (14)

This solution is not a function of $N$ due to the choice of the initial
condition that $P_{cn}(0)=\delta_{c,0}\delta_{n,1}$.

We observe that $\rho(\infty)$ in Eq. (14) is a decreasing function of $p$ but
an increasing function of $q$ [Fig. 1(c)], already partly implying the
qualitatively similar behavior to the simulation results of the original March
model, i.e., Fig. 1 in Ref. [1]. The larger $q$ leads to the more correct
belief of agents about the reality, which is easily understood by considering
that $q$ is the learning probability of the code from superior agents. On the
other hand, the effect of $p$ on $\rho(\infty)$ is not straightforward to
understand. It is because the large value of $p$ speeds up not only the
probability flow to the state $(0,0)$ but also to the state $(1,N)$. It means
that the large $p$ always helps spread the code’s belief to agents whether the
code’s belief is correct or not. When the code’s belief is incorrect, the
large $p$ increases the amount of flow to the state $(0,0)$ [Fig. 1(a)]. In
contrast, when the code’s belief is correct, the large $p$ does not increase
the amount of flow to the state $(1,N)$, but only speeds up the flow. As we
focus on the steady behavior, such an asymmetric role of $p$ leads to the
decreasing behavior of $\rho(\infty)$ as a function of $p$. Such a behavior
has been interpreted that slow socialization allows for longer exploration,
resulting in better organizational learning [1].

For general initial conditions, one can estimate $\rho(T)$ for a sufficiently
large $T$ by iterating the master equation in Eq. (7) for a given initial
condition $\\{P_{cn}(0)\\}$. For a demonstration, we consider a system of
$N=40$ agents and the initial condition that
$P_{cn}(0)=\frac{1}{2}\delta_{n,N/2}$ for each $c\in\\{0,1\\}$. We estimate
the value of $\rho(T)$ at the first time step $T$ when
$|\rho(T)-\rho(T-1)|<10^{-4}$ is satisfied. From the results shown in Fig.
1(d), we find that $\rho(T)$ is a decreasing function of $p$ but an increasing
function of $q$, showing the same tendency as the solution of $\rho(\infty)$
in Eq. (14) for the simpler initial condition.

These exact and numerical results are supported by agent-based simulations. We
perform the simulations of the model using the rules in Eqs. (1) and (2)
together with Eq. (3) for the system with $N=40$ agents with the mentioned
initial conditions. Firstly, the initial condition with
$P_{cn}(0)=\delta_{c,0}\delta_{n,1}$ used for the analysis is realized in the
simulation such that only one agent has an initial belief of $1$, while all
other agents and the code have the belief of $0$. Secondly, as for the initial
condition with $P_{cn}(0)=\frac{1}{2}\delta_{n,N/2}$ for each $c\in\\{0,1\\}$,
we set $\sigma_{i}(0)=1$ for $i=1,\ldots,20$ and $\sigma_{i}(0)=0$ for the
rest of agents, while the value of $c(0)$ is randomly chosen from $\\{0,1\\}$
with equal probabilities. Eventually every run ends up with one of absorbing
states, implying that $n(\infty)=0,N$. For each pair of $p$ and $q$, we take
the average of $n(\infty)/N$ over $10^{4}$ different runs to get the value of
$\rho(\infty)$ in Eq. (5). Such averages are shown with symbols in Fig. 1(c,
d), which are indeed in good agreement with analytical and numerical
solutions, respectively.

###  III.2 Heterogeneous learning in a closed system

Next, we study the heterogeneous version of the model in a closed system with
$r=1$ by using two distinct values of learning probability, i.e., by setting

$\displaystyle\begin{split}&p^{(i)}=p_{1}\ \textrm{for}\ i=1,\ldots,N_{1},\\\
&p^{(i)}=p_{2}\ \textrm{for}\ i=N_{1}+1,\ldots,N,\end{split}$ (15)

where $1\leq N_{1}\leq N-1$ [Fig. 2(a)]. Agents with the larger (smaller)
learning probability among $p_{1}$ and $p_{2}$ can be called fast (slow)
learners [1]. The state of the system at each time step $t$ can be summarized
in terms of the code’s belief $c(t)$, the number of agents with $p_{1}$ whose
belief is $1$, which we denote by $n(t)\in\\{0,\ldots,N_{1}\\}$, and the
number of agents with $p_{2}$ whose belief is $1$, which we denote by
$m(t)\in\\{0,\ldots,N_{2}\\}$. Here $N_{2}\equiv N-N_{1}$. Then the expected
density of agents with belief of $1$ is given as

$\displaystyle\rho(t)\equiv\left\langle\frac{n(t)+m(t)}{N}\right\rangle,$ (16)

which can also be interpreted as the expected belief of a randomly chosen
agent.

Similarly to the homogeneous version of the model, the master equation reads

$\displaystyle
P_{cnm}(t+1)=\sum_{c^{\prime}n^{\prime}m^{\prime}}W_{c^{\prime}n^{\prime}m^{\prime}\to
cnm}P_{c^{\prime}n^{\prime}m^{\prime}}(t),$ (17)

where the transition probabilities are written as

$\displaystyle\begin{split}&W_{000\to
cnm}=\delta_{0,c}\delta_{0,n}\delta_{0,m},\\\ &W_{0n^{\prime}m^{\prime}(\neq
00)\to 0nm}=\begin{cases}{n^{\prime}\choose
n}p_{1}^{n^{\prime}-n}\bar{p}_{1}^{n}{m^{\prime}\choose
m}p_{2}^{m^{\prime}-m}\bar{p}_{2}^{m}\bar{q}&\textrm{if}\ n^{\prime}\geq n\
\&\ m^{\prime}\geq m,\\\ 0&\textrm{otherwise},\end{cases}\\\
&W_{0n^{\prime}m^{\prime}(\neq 00)\to 1nm}=\begin{cases}{n^{\prime}\choose
n}p_{1}^{n^{\prime}-n}\bar{p}_{1}^{n}{m^{\prime}\choose
m}p_{2}^{m^{\prime}-m}\bar{p}_{2}^{m}q&\textrm{if}\ n^{\prime}\geq n\ \&\
m^{\prime}\geq m,\\\ 0&\textrm{otherwise},\end{cases}\\\
&W_{1n^{\prime}m^{\prime}\to 0nm}=0\ \forall n^{\prime},m^{\prime},n,m,\\\
&W_{1n^{\prime}m^{\prime}\to 1nm}=\begin{cases}{N_{1}-n^{\prime}\choose
N_{1}-n}p_{1}^{n-n^{\prime}}\bar{p}_{1}^{N_{1}-n}{N_{2}-m^{\prime}\choose
N_{2}-m}p_{2}^{m-m^{\prime}}\bar{p}_{2}^{N_{2}-m}&\textrm{if}\ n^{\prime}\leq
n\ \&\ m^{\prime}\leq m,\\\ 0&\textrm{otherwise}.\end{cases}\end{split}$ (18)

Here we have used $\bar{p}_{1}\equiv 1-p_{1}$ and $\bar{p}_{2}\equiv 1-p_{2}$.
It is obvious that there are two absorbing states, i.e., $(0,0,0)$ and
$(1,N_{1},N_{2})$, implying that $P_{000}(\infty)+P_{1N_{1}N_{2}}(\infty)=1$.
Calculating Eq. (17) recursively with any initial condition
$\\{P_{cnm}(0)\\}$, one can in principle obtain $P_{cnm}(t)$ for any $c$, $n$,
$m$, and $t$, hence $\rho(t)$ in Eq. (16).

Figure 2: (a) Schematic diagram of the heterogeneous learning model in a
closed system with $r=1$, showing how the code and agents learn from each
other with learning probabilities $p_{1}$, $p_{2}$, and $q$. (b) Transition
structure between states of the system. Each state is denoted by $(c,n,m)$,
where $c$ is the belief of the code and $n$ ($m$) is the number of agents with
learning probability $p_{1}$ ($p_{2}$) whose belief is $1$. Self-loops and
states with $c=1$ are not shown. (c) Analytic solutions of $\rho(\infty)$ in
Eq. (22) (lines) for the case with $N=40$ ($N_{1}=N_{2}=20$) and an initial
condition that $P_{cnm}(0)=\delta_{c,0}\delta_{n,1}\delta_{m,1}$. (d)
Numerical estimation of $\rho(T)$ using the master equation in Eq. (17)
(lines) for the case with $N=40$ ($N_{1}=N_{2}=20$) and an initial condition
that $P_{cnm}(0)=\frac{1}{2}\delta_{n,N_{1}/2}\delta_{m,N_{2}/2}$ for each
$c\in\\{0,1\\}$. Here $T$ is the first time step satisfying
$|\rho(T)-\rho(T-1)|<10^{-4}$. (e) Numerical estimation of expected beliefs of
the code (“code” in the figure), fast-learning agents with $p_{1}=0.9$
(“fast”), slow-learning agents with $p_{2}=0.1$ (“slow”), and all agents
(“avg”) at $t=3$ using the master equation in Eq. (17) for the case with
$N=40$ and $N_{1}=2,4,6,\ldots,38$ (lines). We use an initial condition that
$P_{cnm}(0)=\frac{1}{2}\delta_{n,N_{1}/2}\delta_{m,N_{2}/2}$ for each
$c\in\\{0,1\\}$ and for each $N_{1}$. In panels (c–e), simulation results are
shown with symbols, each symbol was averaged over $2\times 10^{4}$ different
runs, and standard errors are omitted as they are smaller than symbols.

As the simplest yet nontrivial case, let us consider an initial condition that
$P_{cnm}(0)=\delta_{c,0}\delta_{n,1}\delta_{m,1}$. Denoting

$\displaystyle\begin{split}&\alpha_{i_{1}i_{2}}\equiv W_{0i_{1}i_{2}\to
0i_{1}i_{2}},\\\ &\beta_{i_{1}i_{2},j_{1}j_{2}}\equiv W_{0i_{1}i_{2}\to
0j_{1}j_{2}}\ [(i_{1},i_{2})\neq(j_{1},j_{2})],\end{split}$ (19)

the master equations for $P_{011}$, $P_{010}$, $P_{001}$, and $P_{000}$ [Eq.
(17)] are written as follows:

$\displaystyle P_{011}(t+1)=$ $\displaystyle\alpha_{11}P_{011}(t),$
$\displaystyle P_{010}(t+1)=$
$\displaystyle\alpha_{10}P_{010}(t)+\beta_{11,10}P_{011}(t),$ $\displaystyle
P_{001}(t+1)=$ $\displaystyle\alpha_{01}P_{001}(t)+\beta_{11,01}P_{011}(t),$
(20) $\displaystyle P_{000}(t+1)=$ $\displaystyle
P_{000}(t)+\beta_{11,00}P_{011}(t)+\beta_{10,00}P_{010}(t)$
$\displaystyle+\beta_{01,00}P_{001}(t),$

where $\alpha_{11}=\bar{p}_{1}\bar{p}_{2}\bar{q}$,
$\alpha_{10}=\bar{p}_{1}\bar{q}$, $\alpha_{01}=\bar{p}_{2}\bar{q}$,
$\beta_{11,10}=\bar{p}_{1}p_{2}\bar{q}$,
$\beta_{11,01}=p_{1}\bar{p}_{2}\bar{q}$, $\beta_{11,00}=p_{1}p_{2}\bar{q}$,
$\beta_{10,00}=p_{1}\bar{q}$, and $\beta_{01,00}=p_{2}\bar{q}$. After some
algebra, one obtains

$\displaystyle P_{000}(\infty)=$
$\displaystyle\frac{\beta_{11,00}}{1-\alpha_{11}}+\frac{\beta_{11,10}\beta_{10,00}}{(1-\alpha_{11})(1-\alpha_{10})}$
$\displaystyle+\frac{\beta_{11,01}\beta_{01,00}}{(1-\alpha_{11})(1-\alpha_{01})},$
(21)

leading to

$\displaystyle\rho(\infty)=1-\frac{p_{1}p_{2}\bar{q}}{1-\bar{p}_{1}\bar{p}_{2}\bar{q}}\left[1+\frac{\bar{p}_{1}\bar{q}}{1-\bar{p}_{1}\bar{q}}+\frac{\bar{p}_{2}\bar{q}}{1-\bar{p}_{2}\bar{q}}\right].$
(22)

This result is not a function of $N_{1}$ and $N_{2}$ due to the choice of the
initial condition that $P_{cnm}(0)=\delta_{c,0}\delta_{n,1}\delta_{m,1}$. It
is straightforward to prove that setting $p_{1}=p_{2}=p$ reduces the solution
in Eq. (22) to the solution of the homogeneous model with the initial
condition that $P_{cn}(0)=\delta_{c,0}\delta_{n,2}$.

To demonstrate the effect of heterogeneous learning on $\rho(\infty)$ in Eq.
(22), we parameterize $p_{1}=p+\delta$ and $p_{2}=p-\delta$ with non-negative
$\delta$ and $p\in(\delta,1-\delta]$. Here $\delta$ controls the degree of
heterogeneity of agents. As shown in Fig. 2(c), the larger $\delta$ leads to
the higher values of $\rho(\infty)$ in Eq. (22) for the entire range of $p$,
which is consistent with the simulation results of the original March model,
i.e., Fig. 2 in Ref. [1]. Such behaviors can be essentially understood by
comparing the transition probability $W_{011\to 000}=p_{1}p_{2}\bar{q}$ in the
heterogeneous model to its counterpart $W_{02\to 00}=p^{2}\bar{q}$ in the
homogeneous model [Eq. (8)] to get

$\displaystyle\frac{W_{011\to 000}}{W_{02\to
00}}=1-\frac{\delta^{2}}{p^{2}}\leq 1.$ (23)

It implies that for positive $\delta$ the probability flow to the absorbing
state $(0,0,0)$ in the heterogeneous model is always smaller than the flow to
the absorbing state $(0,0)$ in the homogeneous model, hence the larger
$\rho(\infty)$ for the heterogeneous model than for the homogeneous model. We
also remark that the ratio in Eq. (23) gets closer to $1$ for the larger value
of $p$, hence the smaller gap between the heterogeneous and homogeneous
models. Such an expectation is indeed the case as depicted in Fig. 2(c).

For general initial conditions, we numerically estimate $\rho(T)$ for a
sufficiently large $T$ by iterating the master equation in Eq. (17) for a
given initial condition $\\{P_{cnm}(0)\\}$. For a demonstration, we consider a
system of $N=40$ agents ($N_{1}=N_{2}=20$) and the initial condition that
$P_{cnm}(0)=\frac{1}{2}\delta_{n,N_{1}/2}\delta_{m,N_{2}/2}$ for each
$c\in\\{0,1\\}$. We estimate the value of $\rho(T)$ at the first time step $T$
when $|\rho(T)-\rho(T-1)|<10^{-4}$ is satisfied. From the results shown in
Fig. 2(d), we find that $\rho(T)$ has higher values for more heterogeneous
systems.

We also perform the agent-based simulations of the heterogeneous model for the
system with $N=40$ agents ($N_{1}=N_{2}=20$) with the mentioned initial
conditions. Firstly, the initial condition with
$P_{cnm}(0)=\delta_{c,0}\delta_{n,1}\delta_{m,1}$ is realized in the
simulation such that one agent with $p_{1}$ and one agent with $p_{2}$ have an
initial belief of $1$, while all other agents as well as the code have the
belief of $0$. Secondly, as for the initial condition with
$P_{cnm}(0)=\frac{1}{2}\delta_{n,N_{1}/2}\delta_{m,N_{2}/2}$ for each
$c\in\\{0,1\\}$, we set $\sigma_{i}(0)=1$ for $i=1,\ldots,10,21,\ldots,30$ and
$\sigma_{i}(0)=0$ for the rest of agents, while the value of $c(0)$ is
randomly chosen from $\\{0,1\\}$ with equal probabilities. Eventually every
run ends up with one of absorbing states, implying that
$n(\infty)+m(\infty)=0,N$. For each combination of $p_{1}$, $p_{2}$, and $q$,
we take the average of $[n(\infty)+m(\infty)]/N$ over $2\times 10^{4}$
different runs to get the value of $\rho(\infty)$ in Eq. (16). Such averages
are shown with symbols in Fig. 2(c, d), which are indeed in good agreement
with analytical and numerical solutions, respectively.

Finally, we note that our setup for heterogeneous agents is different from
that in the original March model [1]. In the original paper, the heterogeneity
was controlled by the number of agents having $p_{1}$, i.e., $N_{1}$, while
learning probabilities were fixed to be $p_{1}=0.9$ and $p_{2}=0.1$. We test
such original setup using our simplified model both by estimating $\rho(3)$
from the master equations in Eq. (17) and by performing agent-based
simulations up to $t=3$. For a system of $N=40$ agents, we consider
$N_{1}=2,4,6,\ldots,38$ with the initial condition that
$P_{cnm}(0)=\frac{1}{2}\delta_{n,N_{1}/2}\delta_{m,N_{2}/2}$ for each
$c\in\\{0,1\\}$. In addition to the expected belief of all agents in Eq. (16),
we measure the expected belief of fast-learning agents with $p_{1}$, that of
slow-learning agents with $p_{2}$, and that of the code. Results from the
numerical estimation of master equations and from agent-based simulations are
in good agreement with each other as depicted in Fig. 2(e). These results show
the qualitatively same behaviors as in the original March model, i.e., Fig. 3
in Ref. [1].

###  III.3 Homogeneous learning in an open system

Figure 3: (a) Different behaviors of the analytic solution of $\rho(2)$ in Eq.
(28) as a function of $u$ depicted in the plane of $(p,q)$. $\rho(2)$ can have
either a maximum value for $u\in(0,1)$ (red shade, denoted by “max”) or a
minimum value for $u\in(0,1)$ (tin shade, denoted by “min”), or it can
monotonically increase (“increasing”) or decrease (“decreasing”). Four empty
symbols are chosen to demonstrate their different functional forms of
$\rho(2)$ in panel (b). (b) Analytic solutions of $\rho(2)$ in Eq. (28) as a
function of $u$ for several combinations of $p$ and $q$ (lines) with
corresponding simulation results (symbols). We use an initial condition that
$P_{c\sigma}(0)=\delta_{c,0}\delta_{\sigma,1}$. (c) Heatmaps of the analytic
solution of $\rho(\infty)$ in Eq. (34) as a function of $u$ and $v$ for the
case with $p=q=0.5$. Two empty symbols are chosen to demonstrate different
behaviors of $\rho(t)$ in panel (d). (d) Numerical estimation of $\rho(t)$
using the master equation (solid lines) and simulation results (dotted lines)
for $u=0$ and $0.1$ when $p=0.5$, $q=0.5$, and $v=0.02$. In panels (b, d),
each symbol was averaged over $2\times 10^{5}$ different runs. Standard errors
are omitted as they are smaller than symbols.

We finally study the effects of turnover of agents and turbulence of the
external reality on organizational learning. For this, we focus on the
simplest yet nontrivial case with $N=1$, indicating that there is only one
agent in the system. This agent’s belief is denoted by $\sigma(t)$. The case
with general $N>1$ can be studied too within our framework.

We first consider the system with turnover of agents only, while $r(t)=1$ for
all $t$, namely, $u>0$ and $v=0$. Let us denote by $P_{c\sigma}(t)$ the
probability that at time step $t$ the code’s belief is $c$, and the agent’s
belief is $\sigma$. These probabilities satisfy the normalization condition as

$\displaystyle\sum_{c,\sigma\in\\{0,1\\}}P_{c\sigma}(t)=1.$ (24)

They evolve according to the following master equation in discrete time:

$\displaystyle
P_{c\sigma}(t+1)=\sum_{c^{\prime}\sigma^{\prime}}W_{c^{\prime}\sigma^{\prime}\to
c\sigma}P_{c^{\prime}\sigma^{\prime}}(t),$ (25)

where the transition probabilities read using $u^{\prime}\equiv u/2$ and
$\bar{u}^{\prime}\equiv 1-u^{\prime}$ as follows:

$\displaystyle\begin{split}&W_{00\to 00}=\bar{u}^{\prime},\ W_{00\to
01}=u^{\prime},\\\ &W_{01\to
00}=(p\bar{u}^{\prime}+\bar{p}u^{\prime})\bar{q},\ W_{01\to
01}=(\bar{p}\bar{u}^{\prime}+pu^{\prime})\bar{q},\\\ &W_{01\to
10}=(p\bar{u}^{\prime}+\bar{p}u^{\prime})q,\ W_{01\to
11}=(\bar{p}\bar{u}^{\prime}+pu^{\prime})q,\\\ &W_{10\to
10}=\bar{p}\bar{u}^{\prime}+pu^{\prime},\ W_{10\to
11}=p\bar{u}^{\prime}+\bar{p}u^{\prime},\\\ &W_{11\to 10}=u^{\prime},\
W_{11\to 11}=\bar{u}^{\prime}.\end{split}$ (26)

and all other transition probabilities are zero. Note that due to $u>0$, both
$(c,\sigma)=(0,0)$ and $(1,1)$ are no longer absorbing states.

For the steady state, we derive the analytical solution of $\rho(\infty)$ as

$\displaystyle\rho(\infty)=\sum_{c\in\\{0,1\\}}P_{c1}(\infty)=\frac{p+(\frac{1}{2}-p)u}{p+(1-p)u}.$
(27)

This solution is independent of the initial condition, and it is not a
function of $q$ because the change of the code’s belief from $0$ to $1$ is
irreversible; as long as $q>0$, all initial probabilities end up with states
with $c=1$. We also find that $\lim_{u\to 0}\rho(\infty)=1$ and
$\rho(\infty)=1/2$ for $u=1$, both of which are irrespective of $p$. That is,
$\rho(\infty)$ is a decreasing function of $u$ for $u>0$, whereas for $u=0$ it
can have the finite value less than one, e.g., given in Eq. (14), as long as
$p>0$. Thus one can conclude that $\rho(\infty)$ shows an “increasing” and
then decreasing behavior in the range of $0\leq u\leq 1$. This argument is
important to discuss about the optimal turnover of agents that maximizes the
effectiveness of organizational learning.

Next, we focus on the transient dynamics instead of the steady state. Starting
from the initial condition that
$P_{c\sigma}(0)=\delta_{c,0}\delta_{\sigma,1}$, we obtain, e.g., at $t=2$

$\displaystyle\rho(2)=$
$\displaystyle(\tfrac{1}{2}-p)(1-p)u^{2}-(2p^{2}+pq-\tfrac{7}{2}p+1)u$
$\displaystyle+(1-p)^{2}+pq.$ (28)

It turns out that $\rho(2)$ is a quadratic function of $u$, meaning that it
can have either a maximum or minimum value for the range of $u\in(0,1)$, or it
may be a monotonically increasing or decreasing function of $u$, depending on
the choice of $p$ and $q$. Figure 3(a) summarizes such behaviors in the plane
of $(p,q)$, and Fig. 3(b) depicts $\rho(2)$ as a function of $u$ for several
cases of $p$ and $q$. For example, for sufficiently large $q$, $\rho(2)$ is a
monotonically decreasing function of $u$ irrespective of $p$. It implies that
if the code learns fast from the superior agent, the maximal organizational
learning is achieved when there is no turnover. This result can be understood
by considering the fact that the turnover introduces randomness or new
information from outside of the system. In contrast, if the code learns slowly
from the superior agent but the agent learns fast from the code, the maximal
organizational learning is achieved for the largest turnover. In such case,
without turnover, both the code and agent are likely to be stuck in a
suboptimal situation. Thus, the strong turnover may help the system to evade
it. Precisely, we find the increasing and then decreasing behavior of
$\rho(2)$ for $p=0.7$ and $q=0.6$, and the monotonically decreasing behavior
of $\rho(2)$ for $p=0.4$ and $q=0.6$, which are consistent with the results of
the original March model, e.g., Fig. 4 in Ref. [1].

We now consider the effect of turbulence on the organizational learning in the
presence of turnover of agents. For this, we define an extended system
consisting of both the system and the reality, whose states can be denoted by
$(r,c,\sigma)\in\\{0,1\\}^{3}$. Let us denote by $P_{rc\sigma}(t)$ the
probability that at time step $t$ the reality is $r$, the code’s belief is
$c$, and the agent’s belief is $\sigma$. These probabilities satisfy the
normalization condition as

$\displaystyle\sum_{r,c,\sigma\in\\{0,1\\}}P_{rc\sigma}(t)=1.$ (29)

They evolve according to the following master equation in discrete time:

$\displaystyle
P_{rc\sigma}(t+1)=\sum_{r^{\prime}c^{\prime}\sigma^{\prime}}W_{r^{\prime}c^{\prime}\sigma^{\prime}\to
rc\sigma}P_{r^{\prime}c^{\prime}\sigma^{\prime}}(t).$ (30)

Denoting $v^{\prime}\equiv v/2$ and $\bar{v}^{\prime}\equiv 1-v^{\prime}$ and
using Eq. (26), we get the transition probabilities for Eq. (30) as follows:

$\displaystyle\begin{split}&W_{0c^{\prime}\sigma^{\prime}\to
0c\sigma}=\overline{W}_{c^{\prime}\sigma^{\prime}\to
c\sigma}\bar{v}^{\prime},\\\ &W_{0c^{\prime}\sigma^{\prime}\to
1c\sigma}=\overline{W}_{c^{\prime}\sigma^{\prime}\to c\sigma}v^{\prime},\\\
&W_{1c^{\prime}\sigma^{\prime}\to 0c\sigma}=W_{c^{\prime}\sigma^{\prime}\to
c\sigma}v^{\prime},\\\ &W_{1c^{\prime}\sigma^{\prime}\to
1c\sigma}=W_{c^{\prime}\sigma^{\prime}\to
c\sigma}\bar{v}^{\prime},\end{split}$ (31)

where we have used

$\displaystyle\overline{W}_{c^{\prime}\sigma^{\prime}\to c\sigma}\equiv
W_{1-c^{\prime},1-\sigma^{\prime}\to 1-c,1-\sigma}.$ (32)

As $r(t)$ is no longer constant, the average individual knowledge is obtained
as

$\displaystyle\rho(t)=\sum_{r,c,\sigma\in\\{0,1\\}}\delta_{r,\sigma}P_{rc\sigma}(t).$
(33)

After some algebra, we derive an exact solution of $\rho(\infty)$ for the
steady state as follows:

$\displaystyle\rho(\infty)=$
$\displaystyle[(p^{2}q-p^{2}-\tfrac{5}{2}pq+2p+q-1)u^{2}v^{2}-(p^{2}q+p^{2}-\tfrac{7}{2}pq+2p+\tfrac{3}{2}q-1)u^{2}v-(p-\tfrac{1}{2})qu^{2}-(2p^{2}q-2p^{2}-\tfrac{9}{2}pq+3p$
$\displaystyle+\tfrac{3}{2}q-1)uv^{2}+(2p^{2}q-2p^{2}-\tfrac{9}{2}pq+2p+q)uv+pqu+(p^{2}q-p^{2}-2pq+p+\tfrac{1}{2}q)v^{2}-p(pq-
p-q)v]/[2(p^{2}q-p^{2}$
$\displaystyle-2pq+2p+q-1)u^{2}v^{2}-(2p^{2}q-2p^{2}-5pq+4p+3q-2)u^{2}v+(1-p)qu^{2}-(4p^{2}q-4p^{2}-8pq+6p+3q-2)$
$\displaystyle\times
uv^{2}+(4p^{2}q-4p^{2}-7pq+4p+2q)uv+pqu+(2p^{2}q-2p^{2}-4pq+2p+q)v^{2}-2p(pq-
p-q)v].$ (34)

This analytical solution is depicted as a heatmap in Fig. 3(c) for the case
with $p=q=0.5$. We find that for each value of turbulence $v$ there exists an
optimal turnover probability $u^{*}\in(0,1)$ that maximizes the effectiveness
of organizational learning. Such an optimal turnover probability for a given
$v$ is obtained as

$\displaystyle u^{*}(v)=\frac{v^{2}+3v-\sqrt{2v(v+3)(3v+1)}}{v^{2}-3v-2},$
(35)

which is an increasing function of $v$. It implies that the system is required
to have the larger turnover to adapt to the more turbulent external reality.
Yet the value of $\rho(\infty)$ with the optimal turnover tends to decrease
with $v$ [Fig. 3(c)].

Finally, we look at the transient dynamics of $\rho(t)$ for different values
of $u$ when $p$, $q$, and $v$ are given. We numerically obtain $\rho(t)$ by
iterating the master equation in Eq. (25) using the initial condition that
$P_{rc\sigma(0)}=1/8$ for each state. Numerical results are depicted as solid
lines in Fig. 3(d). The agent-based simulations are also performed using the
initial condition that each of $r(0)$, $c(0)$, and $\sigma(0)$ is randomly and
independently drawn from $\\{0,1\\}$. Simulation results are shown as dotted
lines in Fig. 3(d), which are in good agreement with numerical results. These
results are also qualitatively similar to those in the original March model,
e.g., Fig. 5 in Ref. [1].

##  IV Conclusion

In our work, the celebrated organizational learning model proposed by March
[1] has been simplified, enabling us to explicitly write down the master
equation for the dynamics of the model. We have successfully derived exact
solutions for the simplest yet nontrivial cases and numerically estimated
quantities of interest using the master equations, both results are found to
be in good agreement with agent-based simulation results. Our results help
rigorously understand not only the simplified model but also the original
March model to a large extent.

Our theoretical framework for the simplified March model can be applied to the
original March model as well as variants of the March’s model that incorporate
other relevant factors such as forgetting of the beliefs [3, 11] and direct
interaction and communication between agents in the organization [4, 5, 6].
For modeling the interaction structure between agents, various network models
might be deployed [26, 25, 27, 28]. In conclusion, we expect to gain deeper
insights into the organizational learning using our analytical approach.

###### Acknowledgements.

H.-H.J. acknowledges financial support by the National Research Foundation of
Korea (NRF) grant funded by the Korea government (MSIT) (No.
2022R1A2C1007358).

## References

  * March [1991] J. G. March, Exploration and Exploitation in Organizational Learning, Organization Science 2, 71 (1991).
  * Rodan [2005] S. Rodan, Exploration and exploitation revisited: Extending March’s model of mutual learning, Scandinavian Journal of Management 21, 407 (2005).
  * Blaschke and Schoeneborn [2006] S. Blaschke and D. Schoeneborn, The forgotten function of forgetting: Revisiting exploration and exploitation in organizational learning, Soziale Systeme 12, 99 (2006).
  * Miller  _et al._ [2006] K. D. Miller, M. Zhao, and R. J. Calantone, Adding Interpersonal Learning and Tacit Knowledge to March’s Exploration-Exploitation Model, Academy of Management Journal 49, 709 (2006).
  * Kane and Alavi [2007] G. C. Kane and M. Alavi, Information Technology and Organizational Learning: An Investigation of Exploration and Exploitation Processes, Organization Science 18, 796 (2007).
  * Kim and Rhee [2009] T. Kim and M. Rhee, Exploration and exploitation: Internal variety and environmental dynamism, Strategic Organization 7, 11 (2009).
  * Fang  _et al._ [2010] C. Fang, J. Lee, and M. A. Schilling, Balancing Exploration and Exploitation Through Structural Design: The Isolation of Subgroups and Organizational Learning, Organization Science 21, 625 (2010).
  * Sachdeva [2013] M. Sachdeva, Encounter with March’s Organizational Learning Model, Review of Integrative Business and Economics Research 2, 602 (2013).
  * Schilling and Fang [2014] M. A. Schilling and C. Fang, When hubs forget, lie, and play favorites: Interpersonal network structure, information distortion, and organizational learning: When Hubs Forget, Lie, and Play Favorites, Strategic Management Journal 35, 974 (2014).
  * Chanda and Ray [2015] S. S. Chanda and S. Ray, Optimal exploration and exploitation: The managerial intentionality perspective, Computational and Mathematical Organization Theory 21, 247 (2015).
  * Miller and Martignoni [2016] K. D. Miller and D. Martignoni, Organizational learning with forgetting: Reconsidering the exploration–exploitation tradeoff, Strategic Organization 14, 53 (2016).
  * Chanda [2017] S. S. Chanda, Inferring final organizational outcomes from intermediate outcomes of exploration and exploitation: The complexity link, Computational and Mathematical Organization Theory 23, 61 (2017).
  * Chanda  _et al._ [2018] S. S. Chanda, S. Ray, and B. Mckelvey, The Continuum Conception of Exploration and Exploitation: An Update to March’s Theory, M@n@gement 21, 1050 (2018).
  * Chanda and Miller [2019] S. S. Chanda and K. D. Miller, Replicating agent-based models: Revisiting March’s exploration–exploitation study, Strategic Organization 17, 425 (2019).
  * Marín-Idárraga  _et al._ [2022] D. A. Marín-Idárraga, J. M. Hurtado González, and C. Cabello Medina, Factors affecting the effect of exploitation and exploration on performance: A meta-analysis, BRQ Business Research Quarterly 25, 312 (2022).
  * Castellano  _et al._ [2009] C. Castellano, S. Fortunato, and V. Loreto, Statistical physics of social dynamics, Reviews of Modern Physics 81, 591 (2009).
  * Acemoglu and Ozdaglar [2011] D. Acemoglu and A. Ozdaglar, Opinion dynamics and learning in social networks, Dynamic Games and Applications 1, 3 (2011).
  * Sen and Chakrabarti [2014] P. Sen and B. K. Chakrabarti,  _Sociophysics: An Introduction_ (Oxford University Press, Oxford, 2014).
  * Sîrbu  _et al._ [2017] A. Sîrbu, V. Loreto, V. D. P. Servedio, and F. Tria, Opinion dynamics: Models, extensions and external effects, in  _Participatory Sensing, Opinions and Collective Awareness_, edited by V. Loreto, M. Haklay, A. Hotho, V. D. Servedio, G. Stumme, J. Theunis, and F. Tria (Springer International Publishing, Cham, 2017) pp. 363–401.
  * Proskurnikov and Tempo [2017] A. V. Proskurnikov and R. Tempo, A tutorial on modeling and analysis of dynamic social networks. Part I, Annual Reviews in Control 43, 65 (2017).
  * Proskurnikov and Tempo [2018] A. V. Proskurnikov and R. Tempo, A tutorial on modeling and analysis of dynamic social networks. Part II, Annual Reviews in Control 45, 166 (2018).
  * Baronchelli [2018] A. Baronchelli, The emergence of consensus: A primer, Royal Society Open Science 5, 172189 (2018).
  * Anderson and Ye [2019] B. D. O. Anderson and M. Ye, Recent advances in the modelling and analysis of opinion dynamics on influence networks, International Journal of Automation and Computing 16, 129 (2019).
  * Noorazar [2020] H. Noorazar, Recent advances in opinion propagation dynamics: A 2020 survey, The European Physical Journal Plus 135, 521 (2020).
  * Barabási and Pósfai [2016] A.-L. Barabási and M. Pósfai,  _Network Science_ (Cambridge University Press, Cambridge, 2016).
  * Borgatti  _et al._ [2009] S. P. Borgatti, A. Mehra, D. J. Brass, and G. Labianca, Network analysis in the social sciences, Science 323, 892 (2009).
  * Newman [2018] M. E. J. Newman,  _Networks_ , second edition ed. (Oxford University Press, Oxford, United Kingdom ; New York, NY, United States of America, 2018).
  * Menczer  _et al._ [2020] F. Menczer, S. Fortunato, and C. A. Davis,  _A First Course in Network Science_ (Cambridge University Press, Cambridge, 2020).