Influence of technological innovations on industrial production: A motif analysis on the multilayer network

TL;DR

This paper uses motif analysis on multilayer networks constructed from patent and trade data to identify technological synergies that signal local industrial capabilities and competitive advantages.

Contribution

It introduces a novel motif-based analysis method on multilayer networks to detect significant technological combinations linked to industrial production.

Findings

Significant technological synergies are identified beyond noise levels.

Robust evidence of technology combinations leading to competitive advantages.

Method provides insights useful for policy-making in industrial innovation.

Abstract

We study whether specific combinations of technological advancements can signal the presence of local capabilities allowing for a given industrial production. To this end, we generate a multi-layer network using country-level patent and trade data, and perform a motif-based analysis on this network using a statistical validation approach derived from maximum entropy arguments. We show that in many cases the signal far exceeds the noise, providing robust evidence of synergies between different technologies that can lead to a competitive advantage in specific markets. Our results can be highly useful for policy makers to inform industrial and innovation policies.

Tables

Table 1. Table 1: Examples of highly significant Λ Λ \Lambda motifs. The p-value is averaged over all year pairs y 1 = y 2 subscript 𝑦 1 subscript 𝑦 2 y_{1}=y_{2} giving Δ ​ y = 0 Δ 𝑦 0 \Delta y=0 . To make this selection, we picked the most significant pairs (individual links) ( t , p ) 𝑡 𝑝 (t,p) and then choose t ′ superscript 𝑡 ′ t^{\prime} within the region where the average p-value was highest.

p-value	$p$, $t$, $t^{\prime }$
$2\cdot {10}^{-4}$	4701: Wood pulp
	C05B: Lime; Magnesia; Slag; Cements
	C09K: Materials for applications not otherwise provided for
$5\cdot {10}^{-4}$	2605: Mineral products
	C21D: Modifying the physical structure of ferrous metals
	F04F: Pumping of fluid by direct contact of another fluid or by using inertia of fluid to be pumped
$5\cdot {10}^{-4}$	2605: Mineral products
	C21D: Modifying the physical structure of ferrous metals
	F04F: Working metallic powder
$5\cdot {10}^{-4}$	8443: Printing machine
	D02H: Mechanical methods or apparatus in the manufacture of artificial filaments
	G01T Measurement of nuclear or x-radiation
$7\cdot {10}^{-4}$	4703: Chemical wood pulp
	D21F: Decorating textiles
	B27C: Planing, drilling, milling, turning, or universal machines
$8\cdot {10}^{-4}$	2605: Mineral products
	C21D: Modifying the physical structure of ferrous metals
	FF15B: Systems acting by means of fluids in general
$2\cdot {10}^{-3}$	4703: Chemical wood pulp
	D21F: Paper-making machines
	F03D: Wind motors
$2\cdot {10}^{-3}$	4703: Chemical wood pulp
	D21F: Paper-making machines
	D06Q: Decorating textiles
$3\cdot {10}^{-3}$	8519: Sound recording or reproducing apparatus
	G10K: Sound-producing devices
	G01T: Capacitors; rectifiers, detectors, switching devices
$3\cdot {10}^{-3}$	8519: Sound recording or reproducing apparatus
	G10K: Sound-producing devices
	G04f: Time-interval measuring

Equations

$$RCA_{ca}=\dfrac{W_{ca}(y)}{\sum_{a^{\prime}=1}^{N_{a}}W_{ca^{\prime}}(y)}\Bigg{/}\dfrac{\sum_{c^{\prime}=1}^{N_{c}}W_{c^{\prime}a}(y)}{\sum_{c^{\prime}=1}^{N_{c}}\sum_{a^{\prime}=1}^{N_{a}}W_{c^{\prime}a^{\prime}}(y)}\,.$$

$$M_{ca}(y)=\left\{\begin{array}[]{ll}1&\mbox{if }RCA_{ca}(y)\geq 1\\ &\\ 0&\mbox{otherwise}\end{array}\right.$$

$$B_{tp}(y_1,y_2)=\sum _{c=1}^{N_c}\frac{M_{ct}(y_1)}{k_t(y_1)}\frac{M_{cp}(y_2)}{k_c^{(p)}(y_2)}.$$

$${\Lambda }_{t{t}^{\prime }}^p(y_1,y_2)=B_{tp}(y_1,y_2)B_{t^{\prime }p}(y_1,y_2).$$

$$S=-\sum _{\tilde{\mathbf{M}}\in \Omega }P(\tilde{\mathbf{M}})\ln P(\tilde{\mathbf{M}}),$$

$$⟨{\tilde{k}}_a{⟩}_{\Omega }=\sum _{\tilde{\mathbf{M}}\in \Omega }P(\tilde{\mathbf{M}}){\tilde{k}}_a(\tilde{\mathbf{M}})=k_a^{\ast }\forall a\in A,⟨{\tilde{k}}_b{⟩}_{\Omega }=\sum _{\tilde{\mathbf{M}}\in \Omega }P(\tilde{\mathbf{M}}){\tilde{k}}_b(\tilde{\mathbf{M}})=k_b^{\ast }\forall b\in B.$$

$$\begin{array}{rl}0=& \frac{\delta }{\delta P(\tilde{\mathbf{M}})}[S+\omega \left(1-\sum _{\tilde{\mathbf{M}}\in \Omega }P(\tilde{\mathbf{M}})\right)+\\ & +\sum _{a\in A}{\mu }_a\left(k_a^{\ast }-\sum _{\tilde{\mathbf{M}}\in \Omega }P(\tilde{\mathbf{M}}){\tilde{k}}_a(\tilde{\mathbf{M}})\right)+\sum _{b\in B}{\nu }_b\left(k_b^{\ast }-\sum _{\tilde{\mathbf{M}}\in \Omega }P(\tilde{\mathbf{M}}){\tilde{k}}_b(\tilde{\mathbf{M}})\right)].\end{array}$$

$$P(\tilde{\mathbf{M}}\,|\,\{\mu_{a}\},\{\nu_{b}\})=e^{-H(\tilde{\mathbf{M}}\,|\,\{\mu_{a}\},\{\nu_{b}\})}\Big{/}Z(\{\mu_{a}\},\{\nu_{b}\})\,,$$

$$H(\tilde{\mathbf{M}}|\{{\mu }_a\},\{{\nu }_b\})=\sum _{a\in A}{\mu }_a{\tilde{k}}_a(\tilde{\mathbf{M}})+\sum _{b\in B}{\nu }_b{\tilde{k}}_b(\tilde{\mathbf{M}}),$$

$$Z(\{{\mu }_a\},\{{\nu }_b\})=e^{\omega +1}=\sum _{\tilde{\mathbf{M}}\in \Omega }{e}^{-H(\tilde{\mathbf{M}}|\{{\mu }_a\},\{{\nu }_b\})}.$$

$$P(\tilde{\mathbf{M}}|\{{\mu }_a\},\{{\nu }_b\})=\prod _{a\in A}\prod _{b\in B}{\pi }_{ab}^{{\tilde{M}}_{ab}}(1-{\pi }_{ab}{)}^{{\tilde{M}}_{ab}},$$

$${\pi }_{ab}=⟨{\tilde{M}}_{ab}{⟩}_{\Omega }=\sum _{\tilde{\mathbf{M}}\in \Omega }{\tilde{M}}_{ab}P(\tilde{\mathbf{M}}|\{{\mu }_a\},\{{\nu }_b\})=\frac{{\eta }_a{\theta }_b}{1+{\eta }_a{\theta }_b}$$

$$-\frac{\partial }{\partial {\mu }_a}\ln Z(\{{\mu }_a\},\{{\nu }_b\})\equiv ⟨{\tilde{k}}_a{⟩}_{\Omega }=k_a^{\ast }\forall a\in A,-\frac{\partial }{\partial {\nu }_b}\ln Z(\{{\mu }_a\},\{{\nu }_b\})\equiv ⟨{\tilde{k}}_b{⟩}_{\Omega }=k_a^{\ast }\forall b\in B.$$

$$\mathcal{L}(\{{\mu }_a\},\{{\nu }_b\})=\ln P(\mathbf{M}|\{{\mu }_a\},\{{\nu }_b\})=\sum _{a\in A}{k}_a\ln {\eta }_a+\sum _{b\in B}{k}_b\ln {\theta }_b-\sum _{a\in A}\sum _{b\in B}\ln (1+{\eta }_a{\theta }_b),$$

$$\left\{\begin{array}[]{ll}\sum_{b\in B}\dfrac{\eta_{a}\,\theta_{b}}{1+\eta_{a}\,\theta_{b}}=k_{a}\qquad\forall a\in A\\ \sum_{a\in A}\dfrac{\eta_{a}\,\theta_{b}}{1+\eta_{a}\,\theta_{b}}=k_{b}\qquad\forall b\in B\end{array}\right.$$

$${\tilde{B}}_{tp}(y_1,y_2)=\sum _{c=1}^{N_c}\frac{{\tilde{M}}_{ct}(y_1)}{{\tilde{k}}_t(y_1)}\frac{{\tilde{M}}_{cp}(y_2)}{{\tilde{k}}_c^{(p)}(y_2)},$$

$${\tilde{\Lambda }}_{t{t}^{\prime }}^p(y_1,y_2)={\tilde{B}}_{tp}(y_1,y_2){\tilde{B}}_{t^{\prime }p}(y_1,y_2)$$

Full text

I Introduction

Technological innovation is the main driver of modern economic growth freeman1989technology ; romer1990endogeneous . It is therefore not surprising that measuring and predicting the potential impact of technological innovations on export competitiveness has been the central issue of many studies in the last forty years krugman1979model ; soete1980impact ; dosi1990 ; verspagen1992 as well as the focus of a general interest in the field of innovation systems nelson1993national . Several empirical studies that tried to measure such effects of national innovativeness on productivity and trade were also carried out, with mixed results griffith2004mapping ; bronwyn2010 . Overall, such academic effort provided a theoretical framework and empirical stylized facts helpful to understand the aggregate effect of innovation in determining the competitive advantage of countries in different markets. Policy makers are however more interested in identifying the specific technologies that are relevant for specific markets junker2015 ; EC2009 —a task that is much harder to deal with in an organic and objective fashion. Indeed, the scientific effort addressing the impact of specific technologies on specific markets has been limited to ad hoc case studies that are difficult to compare Devaraj2000 ; Sher2005 .

A recent paper of ours multilayer deals with this issue using a multilayer network characterization of the innovation system. Specifically, a three-layered network of innovation activities is derived starting from the three bipartite networks describing the scientific, technological and production activities of countries; the connection of the multilayer network represent the conditional probability that the information produced by an innovation activity (e.g., a technological sector) will be used in another innovation activity (e.g., an industrial product category) after a given time. Grounded on previous fundamental studies of Economic Complexity EFC1 ; EFC2 , this is the first attempt to build a representation of the innovation system as a complex multilayer network.

In this paper we generalize this approach by measuring the potential influence that a pair of activities has on another activity. That is, we go beyond the single link analysis and consider the motifs of the multilayer network Milo2002 ; Saracco2015 . For simplicity we limit our analysis to the relationships between (pairs of) technologies and products, which however represent a crucial aspect of the innovation system—since the interaction between different technologies is often a driver of innovation and progress recombinant . As in multilayer , we carry out a statistical validation of our results against a null network models derived from maximum entropy principles.

II Materials and methods

In order to build the bi-layered network of technologies and products that will be used in our analysis, we start from the following popular databases.

PATSTAT (www.epo.org/searching-for-patents/business/patstat) collects all the patents by different Patent Offices around the world. The basic units of observation in the dataset is a patent family (i.e., the set of patents with common priorities, that is, referred to the same innovation). Each family is related to the countries of origin of the applicants, and to a (set of) technological codes defined by the International Patent Classification (IPC). We define $W_{ct}(y)$ as the number of patent families associated to IPC code $t$ applied by firms located in country $c$ on year $y$.

BACI export data, recorded by UN COMTRADE (https://comtrade.un.org/), collects the import-export flows (quantified in thousands of current US dollars) among countries in the world, related to production as classified using the Harmonized System 2007 of the World Customs Organization. We define $W_{cp}(y)$ as the monetary value of the overall export of country $c$ for product $p$ during year $y$. Note that we use export data as proxies of (competitive) industrial production, as typically done in the Economic Complexity literature.

In this work we consider a time span of data ranging from 1995 to 2012, for which we have a reliable coverage for both patent and export data. For technologies, we use a 4-digits resolution of IPC codes, resulting in a number of technological sectors $N_{t}$ ranging between 629 and 636. For products, we again use a 4-digit resolution of the Harmonized System, resulting in a number of product categories $N_{p}$ ranging between 1140 and 1176. Finally, the number of considered countries $N_{c}$ varies between 66 and 72. The slight variations of these numbers depend on the particular year considered, and are due to geopolitical changes and periodical re-categorization of technologies and products.

Using these basic data, we can define the Revealed Comparative Advantage (RCA) RCA1965 of a country $c$ on an activity $a$ (which is either a technological sector $t$ or and product category $p$) in a given year $y$:

[TABLE]

Thanks to the RCA we can further define for each year $y$ the binary bipartite networks countries-technologies and countries-products. These are respectively represented by the binary biadjacency matrices $\mathbf{M}^{\mathcal{C,T}}(y)$ and $\mathbf{M}^{\mathcal{C,P}}(y)$ whose elements $M_{ct}(y)$ and $M_{cp}(y)$ are:

[TABLE]

(where again $a$ refers to a technological sector $t$ in the case of the bipartite countries-technologies network and to a product category $p$ in the countries-product case).

Once we have the matrices $\mathbf{M}^{\mathcal{C,T}}(y_{1})$ for the year $y_{1}$ and $\mathbf{M}^{\mathcal{C,P}}(y_{2})$ for the year $y_{2}$, in analogy with multilayer , we can construct the assist matrix $\mathbf{B}^{\mathcal{T\to P}}(y_{1},y_{2})$, whose generic element is defined as:

[TABLE]

In the above expression, $k_{t}(y_{1})=\sum_{c=1}^{N_{c}}M_{ct}(y_{1})$ is the number of countries having technology $t$ in their technological portfolio at year $y_{1}$, and $k_{c}^{(p)}(y_{2})=\sum_{p=1}^{N_{p}}M_{cp}(y_{2})$ is the cardinality of the product basket of country $c$ in year $y_{2}$. As explained in multilayer , $B_{tp}(y_{1},y_{2})$ with $y_{1}\leq y_{2}$ gives the conditional probability that a bit of information produced in the technological sector $t$ in year $y_{1}$ arrives (via a random walk on the coupled bipartite network) at the product category $p$ in year $y_{2}$, through one of the countries having $t$ in its technological basket at $y_{1}$ and $p$ in its product basket at $y_{2}$. The elements $B_{tp}(y_{1},y_{2})$ then represent the weighted links of the bi-layered (or bipartite) network technologies-products, with the former at year $y_{1}$ and the latter at year $y_{2}$. These links are extensively studied in multilayer .

Here we move forward and consider the $\Lambda$ motifs of the bipartite technologies-products network:

[TABLE]

$\Lambda_{tt^{\prime}}^{p}(y_{1},y_{2})$ gives the conditional probability that two bits of information originally located on technologies $t$ and $t^{\prime}$ respectively at year $y_{1}$ both reach product $p$ at year $y_{2}$. In other words, this motif quantifies the joint probability for the co-occurrence in a single country of the pair technology $t$ and product $p$ and of the pair technology $t^{\prime}$ and product $p$, where the two events are considered as independent. Note that while this interpretation of the $\Lambda$ motifs cannot be directly related to “impact” or “causality”, it does go beyond a simpler measure of a (time-dependent) correlation. Note also that the name $\Lambda$ motif comes from the fact that, in the bipartite network technologies-products, this quantity gives the weight of the a $\Lambda$ shaped set of two links having different origin in the technology layer and same end in the product layer Saracco2015 . In principle it is possible to generalize this approach by considering higher-order motifs, e.g., by assessing the influence of a wider group of technologies on a single product. For the sake of simplicity and limits of statistical significance, here we focus on the simple $\Lambda$ motif.

After obtaining the empirical values of the $\Lambda$ motifs from data, we statistically validate them using their probability distribution derived from an appropriate null model, which is schematically defined as follows (see Appendix for a thorough presentation). For each year $y$, we build two statistical ensembles of biadjacency matrices $\mathbf{\tilde{M}}^{\mathcal{C,T}}(y)$ and $\mathbf{\tilde{M}}^{\mathcal{C,P}}(y)$ respectively for the bipartite networks countries-technologies and countries-products. These networks are built such to be maximally random, apart from having the ensemble average of node degrees equal to the values observed in the empirical networks. For node degrees, we mean both the technological diversification of countries $\tilde{k}_{c}^{(t)}(y)=\sum_{t}\tilde{M}_{ct}^{y}$ and the technologies ubiquities $\tilde{k}_{t}(y)=\sum_{c}\tilde{M}_{ct}^{y}$ for the countries-technologies network, and both the product diversification of countries $\tilde{k}_{c}^{(p)}(y)=\sum_{p}\tilde{M}_{cp}^{y}$ and the product ubiquities $\tilde{k}_{p}(y)=\sum_{c}\tilde{M}_{cp}^{y}$ for the countries-products network. We choose these quantities as constraints as we want our null model to bear only the information contained in the diversification of countries and the ubiquity of activities, without taking into account the specific pattern of co-occurrences found in the empirical networks. In the spirit of the information theory interpretation of statistical mechanics Jaynes ; Park2004 ; Cimini2018 , the probability measure defining both statistical ensembles of binary bipartite networks is obtained using a constrained entropy maximization approach. The resulting ensembles are known in the literature as Bipartite Configuration Models (BiCM) Saracco2015 . Finally, given the BiCM ensembles for the bipartite networks $\mathbf{\tilde{M}}^{\mathcal{C,T}}(y_{1})$ and $\mathbf{\tilde{M}}^{\mathcal{C,P}}(y_{2})$, we can use Eqs. (3) and (4) appropriately applied to BiCM quantities to derive the probability distribution for the value $\tilde{\Lambda}_{tt^{\prime}}^{p}(y_{1},y_{2})$ in the null model Gualdi2016 ; Saracco2017 . Numerically, we populate the BiCM ensembles by generating $10^{3}$ matrices $\tilde{\mathbf{M}}^{\mathcal{C,T}}(y_{1})$ and $\tilde{\mathbf{M}}^{\mathcal{C,P}}(y_{2})$, and then contract each pair to generate a final ensemble of $10^{3}$ null matrices $\tilde{\mathbf{B}}^{\mathcal{T\to P}}(y_{1},y_{2})$.

In the following we will focus on $\Lambda_{tt^{\prime}}^{p}(\Delta y)$, that is, the average value of $\Lambda_{tt^{\prime}}^{p}(y_{1},y_{2})$ over all pairs of years giving the same difference $y_{2}-y_{1}=\Delta y$. This represents the conditional probability that two bits of information, produced in the same year for a pair of technologies $t$, $t^{\prime}$, reach a product $p$ after $\Delta y$ years. We define the signal $\phi(\Delta y)$ as the fraction of significant $\Lambda_{tt^{\prime}}^{p}(\Delta y)$ (at the $\alpha=0.01$ significance level, according to the probability distribution of the null model) for combinations of $t$, $t^{\prime}$ and $p$ chosen for selected matrix regions. In general we consider a population of motifs equal to 5500 units (see below).

III Results

As first test we report the mean signal $\phi$, that is, the signal averaged over all combinations of $t$, $t^{\prime}$ and $p$. Since the total number of such motifs is extremely large, we choose 5500 motifs at random and take their mean signal as representative of the global average. Figure 1 shows that $\phi$ basically remains within one standard deviation from the noise level $\alpha$, indicating that the mean signal within the data is negligible.

We then report the signal relative to motifs within selected regions of the assist matrix (Fig. 2). Specifically, we choose sub-regions of 11 technologies and 100 products (related to specific technological sectors and product categories), whose total number of $\Lambda$ motifs is 5500 (since $\Lambda_{tt^{\prime}}^{p}=\Lambda_{t^{\prime}t}^{p}$). From Fig. 2 we see that, by selecting coherent sets of technologies and products, the signal is well enhanced: the presence of a pair of technologies in the capability basket of a country can predict whether that country can successfully export a product, and this happens almost independently on the time lag $\Delta y$.

As consistency checks we make two exercises, both reported in Figure 3. Firstly we show that the results we just presented do not depend on the particular resolution used to choose the motifs. Secondly we show that for incoherent technologies and products we indeed get a much lower signal—even lower than the significance level. In this latter case, a significant development of specific technologies corresponds to a low level of export for given products.

We finally provide a few examples of motifs with high signal (i.e., with low p-value). To do that, since the total number of motifs is extremely high so that a complete exploration cannot be performed efficient, we considered the motifs made up of the link pairs $B_{tp}(y_{1},y_{2})$ and $B_{t^{\prime}p}(y_{1},y_{2})$ which are independently those with the highest signal. Table 1 reports some instances of such motifs for the specific choice $\Delta y=0$. The triplets $t,t^{\prime},p$ appearing in the Table seem indeed coherent, and confirm that our method can actually extract meaningful information in an unsupervised way.

IV Conclusions

In this work we provide an effectual way of measuring the combined effect of a set of technologies on one product. In particular we consider lambda motifs, quantifying the paired effect of two technology together. In the process of finding relevant combinations, we highlight several results.

First of all, we show how the combination of multiple technologies has a very different role in different industrial and technological sectors. While technologies related to physics, engineering and chemistry tend to show synergies between them in enhancing the chances of successful export of a product (Figure 2), looking at two generic technologies there is no such effect (Figure 1). This heterogeneity, while expected, is here quantitatively measured. Secondly, we confirm that co-occurrences between technological activities in a country are able to extract information on shared capabilities, to inform policy makers and stakeholders of relevant synergies, like those highlighted in Table 1, for specific export markets.

The mapping provided by our approach for the effects of pairs of technologies on products can represent the fundamental building block for the formulation of a powerful instrument to inform policies and industrial strategies about technology transfer. This operational step will be an important aspect of future research.

V Appendix: the Bipartite Configuration Model (BiCM) and the null model of the Assist matrix

In order to assess the statistical significance of elements of the assist matrices, we resort to a null model for the bipartite matrices $\mathbf{M}^{\cal C,T}(y)$ and $\mathbf{M}^{\cal C,P}(y)$, built by randomly reshuffling their elements (i.e., the network links connecting nodes in the layer ${\cal C}$ of countries to nodes in the layers ${\cal T}$ of technologies and ${\cal P}$ of products, respectively), but preserving the diversification of countries and the ubiquity of the different innovation activities (i.e., the degrees of the nodes in both layers of the bipartite networks countries-technologies and countries-products) . This means that we randomize the signal coming from the network connectivity patters beyond that encoded in the nodes degrees. In order to analytically formulate the null model, avoiding to rely on a conditional uniform graph test Zweig2011 ; Neal2014 , degree constraints are imposed on average, in a way formally similar to what happens for the canonical ensemble in Statistical Mechanics with the constraint on the energy Jaynes . This amounts to set a null hypothesis described by the Bipartite Configuration Model (BiCM) Saracco2015 , which is an extension of the Configuration Model Park2004 to bipartite networks.

In the following, we use symbols with the tilde for quantities assessed on null model configurations, and without the tilde for empirically observed values. From an operational viewpoint, the BiCM null model for any binary biadjacency matrix $\mathbf{M}$, representing a (real) empirical bipartite network with two layers of nodes $a\in A$ and $b\in B$, is built using two main steps:

Through a constrained maximum entropy approach, we define the ensemble $\Omega$ of bipartite networks which are maximally random, apart from the ensemble average of the node degrees on both layers of the bipartite network which are constrained to generic fixed values. Such an ensemble is thus an instance of Exponential Random Binary Graph (ERBG). 2. 2.

In order to determine the ERBG that best represents the empirical bipartite network, we use a maximal likelihood argument showing that the mean values of the node degrees have to be taken equal to the observed ones in the empirical network Squartini2018 : $\langle\tilde{k}_{a}\rangle_{\Omega}=k_{a}$ $\forall a\in A$ and $\langle\tilde{k}_{b}\rangle_{\Omega}=k_{b}$ $\forall b\in B$, where we have indicated with $k$ the observed degrees in the real network and with $\tilde{k}$ the degrees in a generic configuration of the null model. We remind that $k_{a}=\sum_{b}M_{ab}$ and $k_{b}=\sum_{a}M_{ab}$ and analogously for ”tilded” quantities.

Let us start by introducing the ERBG with fixed mean node degrees and let $\tilde{\mathbf{M}}\in\Omega$ be a network configuration in such ensemble and $P(\tilde{\mathbf{M}})$ be its occurrence probability. By implementing the prescriptions from Information Theory and Statistical Mechanics Huang1987 ; Jaynes , the least biased choice of $P(\tilde{\mathbf{M}})$ is the one that maximizes the informational entropy

[TABLE]

subject to the normalization condition $\sum_{\tilde{\mathbf{M}}\in\Omega}P(\tilde{\mathbf{M}})=1$ plus the constraints:

[TABLE]

where $k^{*}_{a}$ $\forall a\in A$ and $k^{*}_{b}$ $\forall b\in B$ are arbitrarily fixed values for the mean degrees of nodes belonging to layers $A$ and $B$. By defining the respective Lagrange multipliers $\omega$, $\{\mu_{a}\}_{a\in{A}}$ and $\{\nu_{b}\}_{b\in{B}}$ (one for each node of the bipartite network), the probability distribution of all configurations $\tilde{\mathbf{M}}\in\Omega$ that maximizes the entropy satisfying at the same time all the constraints, is determined by the following variational equation:

[TABLE]

It is a matter of simple algebra to show that the solution of this equation is:

[TABLE]

where the function $H(\tilde{\mathbf{M}}\,|\,\{\mu_{a}\},\{\nu_{b}\})$ is usually called the Hamiltonian of the graph configurations

[TABLE]

and $Z(\{\mu_{a}\},\{\nu_{b}\})$ is the corresponding partition function

[TABLE]

The above Eqs. (8), (9) and (10) define the network ensemble known as the BiCM model.

Note that, as we have implemented only local constraints, i.e., the mean node degrees, Eq. (8) can be rewritten as the product of single link probability distributions over all pair of nodes belonging respectively to the two different layers Saracco2015 :

[TABLE]

where $\pi_{ab}$ is simply the probability of the link between nodes $a\in A$ and $b\in B$:

[TABLE]

with $\eta_{a}=e^{-\mu_{a}}$ and $\theta_{b}=e^{-\nu_{b}}$. In other words the existence of different links are independent events with respective probabilities which are function only of the Lagrange multipliers associated to the node pairs defining the links. The values of the Lagrange multipliers are determined by the constraints Eqs. (6) which can be rewritten in terms of the derivatives of the partition function:

[TABLE]

Equations (8)-(12) define the generic EBRG with fixed mean degrees of nodes on both layers of the bipartite network.

Once the generic EBRG has been defined, we can move to the step of determining, among all possible EBRGs, the optimal null model for a given real biadjacency matrix (i.e., a bipartite network) $\mathbf{M}$. Equivalently we have to choose the best values for $\{k^{*}_{a}\}_{a\in A}$ and $\{k^{*}_{b}\}_{b\in B}$, i.e., for the Lagrange multipliers $\{\mu_{a}\}_{a\in A}$ and $\{\nu_{b}\}_{b\in B}$, in relation to the connectivity properties of $\mathbf{M}$.

To this aim we write the log-likelihood function Squartini2018

[TABLE]

where $P(\mathbf{M}\,|\,\{\mu_{c}\},\{\nu_{a}\})$ is the probability measure (11) evaluated for a configuration coinciding with the real network $\mathbf{M}$ and $\{k_{a}\}_{a\in A}$ and $\{k_{b}\}_{b\in B}$ are the node degrees of $\mathbf{M}$. The best values for $\{\eta_{a}\}_{a\in A}$ and $\{\theta_{b}\}_{b\in B}$ (or equivalently $\{\mu_{a}\}_{a\in A}$ and $\{\nu_{b}\}_{b\in B}$) are therefore obtained by maximizing such log-likelihood in these parameters. It is simple to show that this amounts to solve the system of $\|A\|+\|B\|$ equations in $\|A\|+\|B\|$ unknowns (where $\|A\|$ and $\|B\|$ simply indicate the number of nodes respectively in the two-layers $A$ and $B$ of $\mathbf{M}$):

[TABLE]

which exactly amounts to choose $k^{*}_{a}=k_{a}$ $\forall a\in A$ and $k^{*}_{b}=k_{b}$ $\forall b\in B$. Finally, the null model for a real bipartite network $\mathbf{M}$ is defined by Eqs. (11) and (12) with the Lagrange multipliers set by Eqs. (15). This recipe can therefore be applied to construct an appropriate null model for all empirical bipartite networks $\{\mathbf{M}^{\cal C,T}(y),\mathbf{M}^{\cal C,P}(y)\}_{y_{min}\leq y\leq y_{max}}$ obtained respectively by the PATSTAT and COMTRADE data.

In order to build a null model for the assist matrices $\mathbf{B}^{\mathcal{T\to P}}(y_{1},y_{2})$, and consequently for the $\Lambda$ motifs defined by Eq. (4), we have now to compose the null models for the bipartite networks $\mathbf{M}^{\cal C,T}(y_{1})$ and $\mathbf{M}^{\cal C,P}(y_{2})$. This is done by contracting the two BiCMs for the matrices $\mathbf{M}^{\cal C,T}(y_{1})$ and $\mathbf{M}^{\cal C,P}(y_{2})$ along the country dimension, as for Eq. (3). We have:

[TABLE]

where $\tilde{k}_{t}(y_{1})=\sum_{c}\tilde{M}_{ct}(y_{1})$ and $\tilde{k}_{c}^{(p)}(y_{2})=\sum_{p}\tilde{M}_{cp}(y_{2})$ are respectively the ubiquity of technology $t$ and the product diversification of country $c$ in the two single configurations for the BiCM null models for the two bipartite networks countries-technologies of year $y_{1}$ and countries-products of year $y_{2}$. In other words, starting from the two BiCM ensembles for $\mathbf{M}^{\cal C,T}(y_{1})$ and $\mathbf{M}^{\cal C,P}(y_{2})$ we build by composition an ensemble of configurations of bipartite networks $\Omega^{\mathcal{T\to P}}(y_{1},y_{2})$ with link weights given by Eq. (16). The probability distributions of elements $\tilde{B}_{tp}(y_{1},y_{2})$, describing the null model, can be in principle obtained using exact techniques Gualdi2016 ; Saracco2017 . However, due to the strong non-Gaussianity of such distributions, we adopt a more practical sampling technique: starting from the BiCMs for $\mathbf{M}^{\cal C,T}(y_{1})$ and $\mathbf{M}^{\cal C,P}(y_{2})$, we use Eqs. (11), (12) and (16) to generate null Assist matrices, and populate the related ensemble $\Omega^{\mathcal{T\to P}}(y_{1},y_{2})$ to estimate the full probability distributions. In a similar way, by using the composition Eq. (4) for the ensemble $\Omega^{\mathcal{T\to P}}(y_{1},y_{2})$:

[TABLE]

and averaging over all pairs of years $y_{1}$ and $y_{2}$ with fixed delay $\Delta y=y_{2}-y_{1}$ to get $\tilde{\Lambda}_{tp}(\Delta y)$, we can construct the null distribution of $\Lambda$ motifs for each triple $t,t^{\prime},p$ and delay $\Delta y$.

The generic observed element $\Lambda_{tt^{\prime}}^{p}(\Delta y)$ is then considered statistically significant depending on the $p$-value that we can infer from its distribution under the null hypothesis. The specific threshold for statistical significance and the size of the generated ensemble vary on the exercises performed (as highlighted in the text). In our case we fixed the statistical significance level at $\alpha=0.01$. It is useful to recall that the two choices, the threshold and the size of the ensemble, are not unrelated: the higher the threshold we want to test, the bigger the sample we require. We consequently extracted for each couple of years $y_{1}$ and $y_{2}$ two ensemble of $1000$ configurations/matrices for the two BiCM null models $\mathbf{\tilde{M}}^{\cal C,T}(y_{1})$ and $\mathbf{\tilde{M}}^{\cal C,P}(y_{2})$, and contracting one to one the configurations in the two ensembles, we get $1000$ values of $\tilde{\Lambda}_{tt^{\prime}}^{p}(y_{1},y_{2})$ to finally determine the statistical significance of the observed value $\Lambda_{tt^{\prime}}^{p}(\Delta y)$ as explained in Sect. II.

A final comment on the issue of multiple hypothesis testing is in order here. Since we do many statistical tests at once (one for each motif of the Assist matrices, in order to determine the true significant elements we should use a correction for the significant threshold (such as Bonferroni or False Discovery Rate). However, the signal $\phi(\Delta y)$ is meant to measure the average outcome of a statistical test over all possible motifs in the network, and as such we do not need to implement any correction for the significance threshold.

Acknowledgement

This work was supported by the Italian PNR Project CRISIS-Lab. G.C. also acknowledge support from the EU project CoeGSS (676547).