Censored Regression for Modelling International Small Arms Trading and its "Forensic" Use for Exploring Unreported Trades

02/25/2019 ∙ by Michael Lebacher, et al. ∙ 0

In this paper we use a censored regression model to analyse data on the international trade of small arms and ammunition (SAA) provided by the Norwegian Initiative on Small Arms Transfers (NISAT). Taking a network based view on the transfers, we do not only rely on exogenous covariates but also estimate endogenous network effects, including a measure for reciprocity as well as the weighted out- and indegree. We rely on a spatial autocorrelation (SAR) model with multiple weight matrices. The likelihood is maximized employing the Monte Carlo Expectation Maximization (MCEM) algorithm. The fitted coefficients reveal strong and stable endogenous network effects. Furthermore, we find evidence for substantial path dependence as well as a close connection between exports of civilian and military small arms. Interestingly, the binary occurrence of major conventional weapons exports is positively related to SAA exports, whereas the volumes have a negative impact. The model is then used in a "forensic" manner to analyse latent network structures and identify countries with higher or lower tendency to export or import than mirrored in the data.



There are no comments yet.


This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

1 Introduction

The Small Arms Survey Update 2018 indicates transfers of small arms in 2015 amounting to 5,7 billion (Holtom and Pavesi, 2018, p. 19) with a major share and highest increases in ammunitions (Holtom and Pavesi, 2018, p. 22). Given the often fatal consequences of the availability of these arms - civilian or military - for intrastate conflict and shootings and for interstate war - the absence of empirical evidence for supplier-recipient networks is surprising. A major reason behind this research gap is the notorious data deficiencies due to non-reporting and illicit trafficking (see Holtom and Pavesi, 2018, p. 29-46). In this research article, we aim to analyse the small arms trading network by combining elements of statistical network data analysis, gravity models and forensic statistical analysis.

Starting with the seminal work of Tinbergen (1962), the gravity equation was quickly established as a valuable tool of empirical trade research. The success of the model stems from its intuitive interpretation as well as its surprisingly strong empirical validity and stability, see for example the meta-analysis of the distance effect by Disdier and Head (2008). It is therefore not surprising that the concept was applied to all kinds of trade relations, including the international exchange of arms and military equipment. An early example for these applications is the work of Bergstrand (1992). Although he doubted the suitability of the model for arms trade because of the strong political considerations in this area, the approach was not abandoned in subsequent research. Akerman and Seim (2014) and Thurner et al. (2018) used elements of the gravity model in order to explain whether Major Conventional Weapons (MCW) are exchanged. Martinez-Zarzoso and Johannsen (2017) rely on the gravity framework of Helpman et al. (2008) in order to investigate the influence of economic and political variables on the extensive and intensive margin of MCW trade. The interplay between oil imports and arms exports is investigated using a gravity model in Bove et al. (2018). While the papers above focus on the exchange of MCW, in this paper we investigate transfers of small arms and ammunition (SAA) provided by the Norwegian Initiative on Small Arms Transfers (NISAT). This data is arguably even better suited for a gravity model since small arms are potentially less sensitive to political decision making and many more trade occurrences are recorded.

We propose a network perspective on international SAA exchanges and conceptualize countries as nodes and transfers between them as directed, valued edges. Although gravity models are a regular part of statistical network analysis (Kolaczyk, 2009), endogenous network effects are rarely incorporated in these models. We do so by connecting the idea of gravity models with the spatial autoregressive (SAR) model adjusted to network data. Especially in sociology, SAR models are regularly used in a network context since the early eighties (Dow et al., 1982; Doreian et al., 1984; Doreian, 1989). They are called network autocorrelation models in this strand of literature. More recently, network autocorrelation models became popular in political science applications, see for example Franzese and Hays (2007), Hays et al. (2010) and Metz and Ingold (2017). Here, it is assumed that actors with certain characteristics are embedded in a network and this embedding leads to contagion and/or spillover effects transmitted through the edges that relate the actors (Leenders, 2002). Hence, one presumes that the characteristics of actors are correlated because a specific social, political or economic mechanism is connecting them. Note that these models are not gravity models in the original sense since the outcome is related to the nodes, and the edges only represent indicators for node dependence. In this paper, we are interested in the dependence among the transfers (instead of the actors) and account for outdegree, indegree, reciprocity and exogenous covariates. A similar model in a non-network context is the spatial gravity model (LeSage and Pace, 2008), that accounts for spatial dependence of the exporter, the importer as well as for the spatial importer-exporter dependence.

Contrary to the typical structure of trade data we observe a high degree of reported non-trade in SAA. In other words, the trade network has a large percentage of zero entries which means that no trade between two states was recorded. To accommodate this zero inflation problem we employ a censored SAR model that can be fitted using the Monte Carlo Expectation Maximization (MCEM) algorithm (Dempster et al., 1977; Wei and Tanner, 1990). There are already several similar EM-based approaches that have been pursued. For instance Suesse and Zammit-Mangion (2017) use the EM algorithm in spatial econometric models, Schumacher et al. (2017) apply an EM-based application to a censored regression model with autoregressive errors, and Vaida and Liu (2009) utilize EM estimation in a censored linear mixed effects model. In Augugliaro et al. (2018) a similar estimation procedure is employed in the context of fitting a graphical LASSO to genetic networks.

While the model application per se provides new insights into SAA trading, another important objective is to make use of the model to explore the validity of reported zero trades. This mirrors a ”forensic” intention in that we want to explore whether unreported trades are likely to have happened based on the fitted model. Though this idea proceeds in the direction of forensic statistics and forensic economics (Aitken and Taroni, 2004; Zitzewitz, 2012) our goal is apparently less ambitious. We do not aim to provide statistical evidence that some states are under-reporting but we do want to explore potential under-reporting by utilizing the fitted network model.

This paper is organized as follows: after presenting the data in Section 2 we explain the model and show how to proceed with estimation and inference in Section 3. In Section 4

the results of the censored regression analysis are given and Section

5 provides the ”forensic” analysis. Section 6 concludes our approach and results.

2 Data Description

Since 2001, the Geneva-based Small Arms Survey specializes on documenting the international flows of the respective products. However, only the Norwegian Initiative on Small Arms Transfers (NISAT, see Marsh, 2017) provides truly relational data necessary for applying network analysis. The NISAT database provides relational information on the trade of small arms, light weapons and ammunition (see also Marsh and McDougal, 2016). This information is collected from different sources as described in Haug et al. (2002). Although NISAT represents the most reliable source of data regarding the exchange of small arms and light weapons, there is nevertheless an enormous amount of uncertainty inherent to arms trading data. This is especially true for light weapons where data quality and availability is partly very poor. Therefore, we restrict our analysis to small arms and the associated ammunition (SAA). See Table 2 in Annex A for the types of small arms and ammunition included in the dataset. Note, that the NISAT database also contains data on sporting guns, which we excluded from the dataset since we are particularly interested in the export of small arms with potential military value. Actually, we will rely on transferred sporting guns volumes later as a useful explanatory variable.

Figure 1: Binary SAA trade network for the most relevant countries in 1992 (left) and 2014 (right). Countries are indicated by grey nodes and directed edges among them indicate arms transfers.

In the remaining dataset, more than 86 000 SAA transfers are recorded for the years 1992-2014, providing the exporting country, the importing country as well as the transferred arms category. The value of the export is measured in constant 2012 USD. In order to make estimation feasible, we restrict our analysis to a subnetwork and select those countries that account for the major share of the SAA trade activity. The resulting 59 countries (see Annex A, Table 1) account for (depending on the year) of the total transfer volume and have participated in arms trade at least once in each year under study. Hence, we investigate the ”core” of the international small arms trade network, balancing the trade-off between the number of countries included, the share of trade volume and the density of the subnetworks. In Figure 1, we show two binary networks for 1992 and 20014, with the countries represented as nodes and the arms transfers as directed edges among them.

Figure 2: Aggregated exports (left) and density (right) in the SAA trade networks of the 59 most relevant countries. Countries with the highest export volume (Germany DEU, Italy ITA and United States USA) are highlighted.

In the left panel of Figure 2 we show the aggregated exports for the most important exporters United States (USA), Germany (DEU) and Italy (ITA) together with the exports of the 56 other countries (other). On the right hand side of Figure 2 we present the density, defined as the sum of existent edges divided by the number of potential edges. Although the network can be without doubt described as a dense one (as compared to the density of social networks), the density is smaller than in the beginning and remains below in the subsequent recent years.

3 Regression Model

3.1 General Model

Let represent a network of transfers at the discrete time points . At each time point consists of nodes and directed, continuous real valued edges, with diagonal elements left undefined. We set

as the row wise vectorization of

, excluding the diagonal elements. In the following, we suppress the time index for ease of the notation and assume (after some suitable transformation) that in each time point follows the autoregressive network model


with being a -dimensional parameter vector for the design matrix . The matrices are row-normalized weight matrices representing linear endogenous network effects, with parameters as their strength. Model (1) is usually known as spatial autoregressive (SAR) model and we refer to LeSage and Pace (2009) for a more detailed discussion and to Lacombe (2004) or LeSage and Pace (2008) for similar models with multiple weight matrices. Standard software implementations that allow for a likelihood based estimation of the model are mostly restricted to the special case with , for example in the R package spdep (Bivand et al., 2013; Bivand and Piras, 2015). The package tnam by Leifeld et al. (2017) allows for multiple weight matrices but is based on pseudo-likelihood estimation and therefore valid only if the weight matrices exclusively apply to exogenous covariates.

Model (1) can be rewritten as

where the dependence on the -dimensional parameter vector is made explicit for notational clarity. Similar as in Besag (1974) and given that all edges in the network are observed, their distribution is given by


The parameter space of this model is restricted such that

is non-singular, which is ensured if the eigenvalues of

are real valued and greater than zero.

3.2 Censored Regression Model

The above model is not directly applicable to our data since a large proportion of the SAA trade values is zero, expressing no (reported) SAA trade between countries and . We therefore adapt model (1) towards a utility model with censored observations. For each potential transfer from to

there exists a utility of the transfer. This utility, however, only materializes in a transfer if it is higher than a certain threshold. Therefore, we assume that the probability model (

2) applies to a network of partly observed latent variables, say . The relation among and is given by for and some threshold. Accordingly, we now set and label the observed utility and the unobserved ones as . A reordering according to the observational pattern of gives

where and . Since the density of the network (see Figure 2) is roughly between and in all years, is always substantially larger than . The mean-covariance structure is given by

In the following, we will denote all reordered matrices in the notation with double subscripts, i.e.  refers to the submatrix of where only interactions of observed variables enter.

3.3 Monte Carlo EM Estimation

In order to estimate the unknown parameter vector , we employ the EM algorithm (Dempster et al., 1977). The complete log-likelihood is simply derived from (2). The respective score vector and Hessian matrix can be found in the Annex B. We are interested in maximizing the observed log-likelihood , where the first part is simply the multivariate normal density of the observed transfers. The second part equals

where and

are the first and second conditional moments. Because

is greater than in each year, the observed log-likelihood is numerically hard to evaluate (and even more so to maximize) with state of the art software implementation. As a solution, we apply the EM algorithm and maximize iteratively. The observation space is given by


3.4 E-Step

The E-Step essentially boils down to calculating the first two moments of a multivariate normally distributed variable


with restriction from (3) applied to . Let those truncated moments be and and define


as the vector that contains the observed values as well as the conditional expectation of the non-observed ones. Given the two moments, we can calculate the conditional expectation of the quadratic form (see Mathai and Provost, 1992):


Then, the function to maximize in the M-step is given by


In order to find the first and second moment of a truncated multivariate normally distributed variable, Vaida and Liu (2009) use the results of Tallis (1961)

on the moment generating function to provide closed form expressions of the E-Step. This, however, is not practicable in our setting as (a) software implementations of a multivariate normal distribution function are overstrained by the high dimension of our problem (the standard package in

R, mvtnorm by Genz et al. (2016) is not able to process dimensions higher than ) and (b) as noted by Schumacher et al. (2017), even if the distribution function could be evaluated, the closed form solution is computationally very expensive which leads to infeasible convergence times in most applications with a high number of non observed values. The same is true for the direct calculation using the moment generating function implemented in R by Wilhelm et al. (2012).

A practicable alternative consists in using the Monte Carlo EM (MCEM) algorithm (Wei and Tanner, 1990) where intractable expectations are replaced by sample based approximations. In our specific case we use the R package TruncatedNormal by Botev (2017) in order to draw from the truncated multivariate normal distribution.

3.4.1 M-Step

It is numerically more efficient to reduce the log-likelihood to a profile log-likelihood by first maximizing with respect to and and then with respect to . Using the derivatives of (7) with respect to and and defining and as the solutions of the score equations as functions of it follows that


where is the hat matrix. With being a constant we can write the profiled function as


The expressions and have derivatives

Now define

which gives


More details on the practical implementation are provided in the Appendix C. Iteration between the E- and the M-step provide the final estimate

. The variance of

can be calculated using Louis (1982) formula with details provided in Appendix D.

4 Application to the Data

4.1 Covariates

Considering model (1) we need to specify the two major components of the model, namely (a) the covariates included in matrix and (b) the network related correlation structure. We start with the first.
Node Specific Variables: Following standard applications (Ward et al., 2013; Head and Mayer, 2014; Egger and Staub, 2016; Thurner et al., 2018), we include the logarithmic real GDP in constant 2010 USD as a measure for the market size of the exporting and importing country. The data is provided by the World Bank (2017)

. For the two years 1993-1994 no reliable GDP data are available for Serbia, Croatia, Estonia, Latvia, Lithuania and Slovenia, we therefore assume that the GDP remained constant in the first three years for this countries. In order to control for the potential influence of intrastate conflicts we include a binary variable that is one if there is an intrastate conflict in the receiving country in the respective year and zero otherwise. The corresponding data is available from the webpage of the Uppsala Conflict Data Program (

UCDP, 2019).

Edge Specific Distance Measures: Because of the strong empirical evidence that distance is a relevant factor in trade (Disdier and Head, 2008), we include the logarithmic distance between capital cities in kilometres (Gleditsch, 2013). In recent applications of the gravity model to arms trade (Akerman and Seim, 2014; Martinez-Zarzoso and Johannsen, 2017; Bove et al., 2018; Thurner et al., 2018) it is argued that political distance measures must also be included in the gravity equation. We therefore also include the absolute difference of the polity IV index (Marshall, 2017

), ranging from 20 (highest ideological distance) to zero (no ideological distance). Additionally, we include a dummy variable for formal alliances between the exporting and importing country, being one if the two countries have a formal alliance. The data is available from

Correlates of War Project (2017) until 2012 and we assume that the alliances stay constant for the years 2013 and 2014.

Edge Specific Trade Measures: We control for lagged logarithmic SAA transfers by smoothing the past observed trade volume using a five-year moving average. In the years with less then five lagged periods available, the moving average is shortened accordingly. Path dependency arises because of diminishing transaction costs, trust relations, security aspects and potentially interoperability, and is a very important determinant in the MCW trade network (Thurner et al., 2018).

Additionally, we include a five year moving average of logarithmic civilian weapon transfers. The intuition behind that is the reasoning that exports of SAA for military usage and civilian usage might be correlated. This is plausible because countries that export massive amounts of civilian arms also have the capabilities to produce military arms. The data is also provided by NISAT (Marsh and McDougal, 2016). Furthermore it seems plausible that there is a connection between the volume of small arms traded and the volume of MCW traded. MCW transfers are recorded by the Stockholm International Peace Research Institute (SIPRI) and measured in so called trend indicator values (TIV). This measure represents the military value and the production costs of the transferred products. For detailed explanation of the data and the TIV see SIPRI (2017b, a) and Holtom et al. (2012). We include a dummy variable that is one if there was an MCW transfer from country to in the actual year or in the four preceding years, zero otherwise. Additionally, we use the logarithmic sum of the exported TIV volumes in the actual year and the four preceding ones.

Figure 3: Schematic representation of linear network effects. The focal edge in dashed grey.

4.2 Network Structure

Next we specify the network specific effects represented by matrices in model (1). We include three effects which are explained subsequently and visualized in Figure 3.

Reciprocity: The reciprocity effect measures whether the export volume from country to country increases in the export volume from to . In the given context it is a plausible assumption that countries tend to specialize in certain types of small arms and/or ammunition and therefore complement each other with their products. Mutual trade is likely to be encouraged by political partnerships and indicates strategic elements, induced by bilateral agreements. The measure is also investigated in the context of commercial trade (e.g. Garlaschelli and Loffredo, 2005; Barigozzi et al., 2010;Ward et al., 2013). In the arms trade literature, reciprocity is investigated empirically by Thurner et al. (2018), with the finding that this is rather unusual in the context of MCW.

Exporter and Importer Effect: The exporter and the importer effect have their analogies in binary networks and can be interpreted as the valued versions of the outdegree and the indegree. The coefficient of the exporter effect measures whether the transfers going out from a certain exporter are correlated. A positive effect indicates the presence of ”super-exporters”. Contrary, the importer effect measures whether the imports of a certain importer are related, with a positive effect indicating ”super-importers”. The degree structure is a crucial feature of the SAA network because a rather small number of countries accounts for the major share of the trade volume, while a small share of (potentially identical) importing countries accounts for a great amount of the import volume.

Before fitting the model we apply the natural logarithm to the data to justify normality. Hence, if the original trade matrices are given by , the elements of are given by if and are not defined if . Furthermore, we define as the lowest strictly positive value in the network at year and set . That is, the threshold is defined such that at a given time point all transfers below the smallest observed log-transformed transfer in that sample are censored. Utility below the threshold implies that no transfer was carried out or was not recorded. Furthermore, we allow for time-varying coefficients by estimating each time-period separately. This relaxes the unrealistic assumption of time-constant effects for more than years and reduces the computational effort. Given these specifications, the final model is now given by

4.3 Results: Coefficients

Figure 4: Time-series of annually estimated regression coefficients. Shaded areas give

2 standard errors. Colouring according to p-values, green:

, yellow: and red: .

In Figure 4 the time series of coefficients are plotted against time for the years 1993-2014. The shaded areas around the coefficients give two standard error bounds and the colouring of the point estimates reflects the respective significance level, the zero-line is depicted by solid black.

Covariates: The exogenous covariates in the first row and in the second row on the right represent the standard gravity variables logarithmic GDP of the exporter and the importer as well as the geographical distance between them (second row, right panel). Overall, the expected results of the gravity equation hold, except for the logarithmic GDP of the exporter that tends to be insignificant and is close to zero in the most recent years. This is an interesting result, because it highlights the fact, that market size is not a prerequisite for producing and exporting internationally competitive SAA. This finding is in stark contrast to the insights on MCW by Thurner et al. (2018).

The strong negative effect of the geographic distance is significant in all years. Again, this is different as compared to the transfers of MCWs where geopolitical strategy disregards distance.

Regarding the political security measures we find that the presence of a conflict in the importing country (second row, left panel) has a mostly positive but seldom significant effect while the coefficient on the distance in political regimes (third row, left panel) is mostly negative but also often insignificant. The coefficient on the dummy variable for formal alliances (third row, right panel) is positive in the beginning but almost permanently insignificant from 2001 on.

The high and consistently significant coefficients of the lagged moving average (fourth row, left panel) illustrates an important feature of the network, namely path inertia. Intensive transfer relationships in the past, strongly increase the export volume in the present. Similarly we find a strong connection between exporting civilian and military arms (fourth row, right panel). Looking at the relation between SAA and MCW trade we find that having traded MCW (fifth row, left panel) in the actual year or the in the four preceding ones has a strong positive effect - at least until the last two years. However, at the same time the effect of the logarithmic sum of the TIV values (fifth row, right panel) has a negative and mostly significant effect. I.e. rather small transfers of MCW tend to coincide with SAA exports while dyads with high amounts of MCW exchange tend to transfer small arms to a relatively lower degree.

Network Structure: On the right panel in the sixth row of Figure 4 the coefficients for reciprocity are shown. The coefficients remain almost constant and positive with values between and . As the coefficients often changes from significant to insignificant we infer that there is at least a tendency that mutuality increases the volume of arms exchanged.

The strongest endogenous effect is the exporter effect (bottom row, left panel) with coefficients that are consistently positive and significant. This indicates that the transfers stemming from the same exporter are indeed highly correlated and reflects the existence of ”super-sellers” like the United States, Germany, Brazil or Italy. On the other hand, we also find a stable, positive and significant importer effect (bottom row, right panel). The fact that the two coefficients on the exporter and the importer effect are much higher than the reciprocity effect provides structural information about heterogeneity in the network. Being a strong exporter or sending to a strong importer increases the export volume more than simply having imported high amounts from the respective partner.

5 ”Forensic” Statistical Analysis of the latent Utility Network

Our model rests on the assumption that the SAA network is determined by a latent utility network

. Based on the joint distribution (

2) we can in fact estimate the probability of being greater than the threshold , given the covariates, the endogenous effects and the rest of the network. In order to do so, let represent the -dimensional vector that contains the realized and the expected values of the latent variables, except the entry that corresponds to the transfer from to . Because we are interested whether some latent transfers could have realized according to the model, we form the expectations without the restriction that the latent transfers must be smaller than . Based on this, we define the conditional probability of a specific latent transfer being greater than the threshold by

See Appendix E for the derivation of . By construction, is high for transfers that are observed in the dataset () and small for most of the transfers that are not within the dataset (). However, we may calculate a high value of , that is a high probability for a realized transfer of arms, despite the data actually indicates . This can be labelled as false positive and we propose to consider this as potential under-reporting. Such a zero-record can happen due to random fluctuation, factors beyond the model as for example historical relationships, or because de-facto existent transfers have not been reported. Vice versa, we may obtain a low value of although . We label this as false negative. Apparently, this requires the fixation of a threshold value for the probabilities. Based on Receiver-Operating Characteristic (ROC) curves, an optimal threshold value can be found using Youden’s statistic (Youden, 1950). This value is optimal in the sense that it allows for a separation such that both, sensitivity and specificity are maximized. This defines the binary network

This network needs is now set into relation with the observed binary SAA trade

Comparing and , we can define the ”forensic” network

which in turn creates two new binary networks

For , the model predicted a transfer that is not present in the dataset, and for , the model did not predict an actual transfer. Therefore, we can label as the false positive network of transfers that are predicted but not realized in the data. Countries with a high out- and/or indegree in the false positive network have more utility from trade relations than mirrored in the data, i.e. it is surprising that the transfers in are not present in the data.

Figure 5: Densities of the false positive network and the false negative network over time.

On the other hand, we refer to as the false negative network

of unpredicted but realized transfers. Here, the model allocates low utility from SAA trade, and the corresponding edge is therefore classified as unlikely, although it is realized.

Figure 6: Ordered box-plot representation of topological network features of the false positive networks for : Eigenvalue centrality (top), outdegree (middle) and indegree (bottom).

In Figure 5 we show the development of the densities of both latent networks.

Figure 7: Ordered box-plot representation of topological network features of the false negative networks for : Eigenvalue centrality (top), outdegree (middle) and indegree (bottom).

The next step is to explore country effects. To do so we evaluate node (i.e. country) specific network topologies of and for each year and summarize the information in box-plots for each country, ordered according to the median of the respective feature. This is shown in Figure 6 for false positives (potential under-reporting) and in Figure 7

for the false negatives. In the first row, we represent the Eigenvector centrality scores. This measure is undirected and constructed such that the centrality of each country is proportional to the sum of the centralities of its trading partners. Hence, countries with high scores have many false positive (false negative) import- and export-relations with many other countries that themselves have many false positive (false negative) import- and export-relations, see e.g.

Csardi and Nepusz (2006). In the middle row, we present the outdegree, that is the number of false positive (false negative) exports for a country. The bottom row in Figures 6 and 7 gives the indegree, that is the number of false positive (false negative) imports. All measures are scaled to take values between and . Countries at the right hand side in the plots of Figure 6 have many false positives and are potentially under-reporting. In Figure 7, the right hand side of the plots mirror many false negatives.

In order to further investigate the measures plotted in Figures 6 and 7 we now focus on the relation between countries to see whether false negatives or false positives occur frequently for particular pairs of countries. To detect persistent patterns in the networks, we check whether false positives or false negatives occur frequently, i.e. counting instances of and for . Denote the aggregated ”forensic” networks as

We look at the distribution of elements of and , which is plotted in Figure 8. On the horizontal axis we show the possible values of the matrix entries, that is the number years where transfers in the forensic networks occur. This ranges from 1 (one false negative or false positive in all years) to 22 (false negatives or false positives in every year from 1993 to 2014). The maximum entry of is thereby less then 22, namely 21, while the maximum value of is 15. On the vertical axis of Figure 8 we show the frequency of the entries of and . Apparently for ”forensic” purposes, large values of are of particular interest, since they report pairs of countries which are likely to under-reporting.

The line in solid black with the ”+” symbols represents and the line in grey with the ”-” symbols represents the false negative network. Additionally, we indicate for both networks the pairs of countries (i.e. sender and receiver) which are of particular interest for ”forensic” purposes. This means for example for an element of that has value 21, that the respective transfer from to is one of the four transfers appeared that appeared times in the false positive network.

Figure 8: Frequency distribution of transfers in the aggregated false positive (, black ”+”) and false negative (, grey ”-”) networks on the vertical axis. Number of years with false positives () or false negatives () on the horizontal axis. Transfers with the most years predicted are indicated in the form ”exporter-importer” in black for and in grey for .

False positive networks : Looking at the Eigenvector centrality scores of Figure 6 on the top provides conclusive results about countries that are central in the network series . Among the countries where arms transfers are potentially under-reported according to the model, are many western European countries as Belgium (BEL), Sweden (SWE), France (FRA), Spain (ESP) and Denmark (DNK). However, the list of presumed under-reporting is headed by Russia (RUS) and Turkey (TUR) but also Brazil (BRA), Israel (ISR) and China (CHN) have high scores. These countries also play a dominant role in Figure 8. In particular, exports from Brazil (BRA) to Russia (RUS), Hungary (HUN), Ukraine (UKR), China (CHN) and Japan (JPN) are likely to be frequently under-reported. Similarly, exports from Russia (RUS) to Belgium (BEL), Cyprus (CYP), Denmark (DNK), Ireland (IRL), Turkey (TUR), Sweden (SWE), Portugal (PRT) and Greece (GRC) are listed. We also find imports of Israel (ISR) from Finland (FIN), Belgium (BEL) and Austria (AUT) as well as exports of Belgium (BEL) to India (IND), Ukraine (UKR), Russia (RUS), Lebanon (LEB), Colombia (COL) and China (CHN).

False negative networks : Among the twelve countries with the highest Eigenvector centrality in Figure 7 is Croatia (HRV) as the only European country. There are however many countries from Asia as Singapore (SGP), India (IND), Thailand (THA), Indonesia (IDN), Malaysia (MYS), South Korea (KOR) and Philippines (PHL). Furthermore, South Africa (ZAF), New Zealand (NZL), Australia (AUS) and Israel (ISR) are among the countries where trade activity is often false negative. For Asian countries as well as for Australia and New Zealand this might mirror the fact that those countries export many SAA to Europe and the United States despite the strongly negative distance effect of the model. Furthermore, this network is very likely to be driven by bilateral agreements and historical developments not covered by the covariates. See for example in Figure 8 the number of false negatives related to the Baltic countries Estonia (EST) Lithuania (LTU) and Latvia (LVA).

It remains to highlight that the constructions of the forensic networks relies on our model with corresponding assumptions and admittedly high degrees of uncertainty. As a consequence, it does not allow for definite statements about actual hidden transfers. However, many of the dyads listed in Figure 8 indeed have either traded massive amounts of civilian arms (e.g. AUT-CHN, BRA-RUS, RUS-BEL, RUS-DNK) or had frequent MCW trade relations (e.g. RUS-CYP) but almost no documented small arms transfers for military usage. Additionally, many of the countries that take central positions in the forensic networks are known for not being very transparent with respect to their SAA exports and imports (see e.g. the small arms transparency barometer).

6 Conclusion

In this paper we have modelled the volumes of international transfers of small arms and ammunition for the years 1992-2014 based on data provided by NISAT. As an analytical tool we combined the gravity model of trade with a modified SAR model that allows to enrich the analysis by endogenous network dependencies, accounting for exporter-related, importer-related and reciprocal dependency among the transfers in the network. Using a censored normal regression model we are able to include information provided by zero-valued transfers. The infeasible likelihood of the censored model is maximized using a Monte Carlo EM algorithm. The fitted model shows strong and stable endogenous network effect, especially related to the sender effect and the receiver effect but also some evidence for reciprocity. Additionally, we find a high coefficient on path dependency and a close connection to the exports of civilian small arms. Conditional on that, the classical gravity hypothesis is confirmed with respect to the GDP of the importer and physical distance but only exceptionally with respect to political distance measures and the GDP of the exporter. This contrasts with the MCW network where distance plays no role, where political similarity and GDP of the exporter have a strong impact (see Thurner et al., 2018). Actually, this difference is plausible, as the technological requirements for the production of small and ammunition are relatively small, and strategic considerations of world-wide acting countries make geographic distances a negligible factor for MCW trade.

Building on our latent utility framework we were able to explore latent utility networks. With the construction of false positive and false negative networks we perform for the first time a forensic approach in this area highlighting especially potentially under-reported exports of Russia and Turkey. We refrain, of course, from making clear-cut assertions. Note that we do not claim to provide unambiguous claims for intentional false reporting. However, we demonstrate that some zero entries in the SAA trading network tend to be not plausible.

Appendix A Descriptives

Country ISO3 Code Country ISO3 Code Country ISO3 Code
Argentina ARG India IND Poland POL
Australia AUS Indonesia IDN Portugal PRT
Austria AUT Ireland IRL Romania ROM
Belgium BEL Israel ISR Russia RUS
Brazil BRA Italy ITA Saudi Arabia SAU
Bulgaria BGR Japan JPN Serbia SRB
Canada CAN Kenya KEN Singapore SGP
Chile CHL South Korea KOR Slovenia SVN
China CHN Kuwait KWT South Africa ZAF
Colombia COL Latvia LVA Spain ESP
Croatia HRV Lebanon LBN Sweden SWE
Cyprus CYP Lithuania LTU Switzerland CHE
Denmark DNK Malaysia MYS Thailand THA
Egypt EGY Mexico MEX Turkey TUR
Estonia EST Netherlands NLD Ukraine UKR
Finland FIN New Zealand NZL Un. Arab Emirates ARE
France FRA Norway NOR United Kingdom GBR
Germany DEU Pakistan PAK United States USA
Greece GRC Peru PER Uruguay URY
Hungary HUN Philippines PHL - -
Table 1: The 59 major exporting and importing countries of the small arms and ammunition dataset with ISO 3 country codes.
Code PRIO Weapons Type Subcategories
200 Small Arms
210 Pistols & Revolvers
230 Rifles/Shotguns (Military)
233 Assault Rifles
234 Carbines
235 Sniper Rifles
237 Semi-automatic Rifles (Military)
239 Shotguns (Military)
240 Machine Guns
243 Sub Machine Guns
245 Light Machine Guns
247 General Purpose Machine Guns
250 Military Weapons
260 Military Firearms
270 Machine Guns All Types
300 Light Weapons
310 Heavy Machine Guns <= 12.7mm
400 Ammunition
415 Small Arms Ammunition
417 Small Calibre Ammunition <= 12.7mm
418 Shotgun Cartridges
Source: nisat.prio.org.
Table 2: Different arms types included in the NISAT dataset with three digit arms category code, weapon type, subcategory and number of transfers in the dataset.

Appendix B Derivatives of the complete log-likelihood

The complete log-likelihood is given by

with score vector


And the corresponding Hessian results in


Where we use Jacobi’s formula (see Magnus and Neudecker, 1988) that allows to express the derivative of a matrix determinant in terms of the derivative of the matrix and its adjugate (). Resulting in

for the third equation in (11). The differentiation of the trace

is used for the sixth and seventh equation in (12).

Appendix C Practical Implementation of the Algorithm

The gradient (10) can be used to maximize (9) by applying the BFGS optimization routine (see Broyden, 1970, Fletcher, 1970, Goldfarb, 1970 and Shanno, 1970). The implementation of the BFGS algorithm in R (R Core Team, 2016) is provided by the base function optim. More computational stability for the maximization of equation (9) is reached by defining as the vector of eigenvalues of and replacing by in equation (9), see Bivand and Piras (2015). The starting value for the algorithm can be found by using a maximum pseudolikelihood estimate (MPLE), using as exogenous covariates in a censored regression model, provided by the R package censReg (Henningsen, 2013). Since the observed log-likelihood cannot be evaluated, we define as the solution of the maximization problem if , otherwise we set and re-iterate until the stopping criteria is satisfied.

Appendix D Approximation of the Fisher Information

Louis (1982) and Oakes (1999) provide formulas for the Fisher information of the observed likelihood. We follow the recommendation of McLachlan and Krishnan (2007), arguing that Louis’s formula is best suited for the MCEM and provides a conservative measure of the standard errors. Therefore, we calculate the observed information based on


Note that the second term of (13) depends not only on the first and second but also on the third and fourth conditional moment of the truncated multivariate normal and cannot be evaluated analytically therefore.

In order to approximate the observed information we are using the results of Robert and Casella (2004, p. 187) and Kang et al. (2013, Section 3.7) that allow for an approximation of the observed information based on the score and Hessian of the complete likelihood. Hence, we can use the results from Annex B for the following procedure.

We draw times potential realizations from the truncated version of (4) using the package TruncatedNormal (Botev, 2017). Those are stored for each draw in a vector , similar as in equation (5).

Then we calculate the score and Hessian from equations (11) and (12) times, where we replace by in each equation and index them by , allowing to calculate the empirical version of (13) by approximating the expectations by means.

This gives an estimator for the observed information. Standard errors are obtained by the square root of the diagonal elements of the inverted approximated matrix.

Appendix E Conditional Probabilities

Based on the fitted coefficients and our model assumptions, we can represent the joint distribution of the latent utility network via a multivariate normal

where and . Given that, define as the -dimensional vector, containing all entries of except . Additionally, for example in the case that is the first entry of , rearrange such that

Then, the conditional distribution of is given by a univariate normal distribution

We are interested in a possible state of the network, where the latent utility is allowed to be greater . Therefore, we insert the expectation for the non-observed utility in and denote this by . Consequently, we can calculate the probability of being greater than using

The probability can be interpreted as the probability that the latent utility of a transfer from country to country is higher than the threshold conditional on the covariates and the remaining network, where no transfer is restricted to be smaller .


We would like to thank Nic Marsh from NISAT for his useful comments and explanations on the data. The project was supported by the European Cooperation in Science and Technology [COST Action CA15109 (COSTNET)]. We also gratefully acknowledge funding provided by the German Research Foundation (DFG) for the project KA 1188/10-1: International Trade of Arms: A Network Approach


  • Aitken and Taroni (2004) Aitken, C. G. and F. Taroni (2004). Statistics and the evaluation of evidence for forensic scientists. Chichester: John Wiley & Sons.
  • Akerman and Seim (2014) Akerman, A. and A. L. Seim (2014). The global arms trade network 1950–2007. Journal of Comparative Economics 42(3), 535–551.
  • Augugliaro et al. (2018) Augugliaro, L., A. Abbruzzo, and V. Vinciotti (2018). L1-penalized censored gaussian graphical model. Biostatistics, https://doi.org/10.1093/biostatistics/kxy043. to appear.
  • Barigozzi et al. (2010) Barigozzi, M., G. Fagiolo, and D. Garlaschelli (2010). Multinetwork of international trade: A commodity-specific analysis. Physical Review E 81(4), 046104.
  • Bergstrand (1992) Bergstrand, J. H. (1992). On modeling the impact of arms reductions on world trade. In C. Isard and W. Anderton (Eds.), Economics of arms reduction and the peace process, pp. 121–142. Amsterdam: Elsevier Science Publishing.
  • Besag (1974) Besag, J. (1974). Spatial interaction and the statistical analysis of lattice systems. Journal of the Royal Statistical Society. Series B 36(2), 192–236.
  • Bivand et al. (2013) Bivand, R., J. Hauke, and T. Kossowski (2013).

    Computing the jacobian in gaussian spatial autoregressive models: An illustrated comparison of available methods.

    Geographical Analysis 45(2), 150–179.
  • Bivand and Piras (2015) Bivand, R. and G. Piras (2015). Comparing implementations of estimation methods for spatial econometrics. Journal of Statistical Software 63(18), 1–36.
  • Botev (2017) Botev, Z. I. (2017). The normal law under linear restrictions: simulation and estimation via minimax tilting. Journal of the Royal Statistical Society: Series B (Statistical Methodology) 79(1), 125–148.
  • Bove et al. (2018) Bove, V., C. Deiana, and R. Nisticò (2018). Global arms trade and oil dependence. The Journal of Law, Economics, and Organization, http://dx.doi.org/10.1093/jleo/ewy007.
  • Broyden (1970) Broyden, C. G. (1970). The convergence of a class of double-rank minimization algorithms 1. general considerations. IMA Journal of Applied Mathematics 6(1), 76–90.
  • Correlates of War Project (2017) Correlates of War Project (2017). Formal interstate alliance dataset, 1648-2012, version 4.1. http://www.correlatesofwar.org/data-sets/formal-alliances. Accessed: 2017-02-06.
  • Csardi and Nepusz (2006) Csardi, G. and T. Nepusz (2006). The igraph software package for complex network research. InterJournal, Complex Systems 1695(5), 1–9.
  • Dempster et al. (1977) Dempster, A. P., N. M. Laird, and D. B. Rubin (1977). Maximum likelihood from incomplete data via the em algorithm. Journal of the Royal Statistical Society. Series B 39(1), 1–38.
  • Disdier and Head (2008) Disdier, A.-C. and K. Head (2008). The puzzling persistence of the distance effect on bilateral trade. The Review of Economics and Statistics 90(1), 37–48.
  • Doreian (1989) Doreian, P. (1989). Models of network effects on social actors. In L. Freeman, D. White, and K. Romney (Eds.), Research Methods in Social Network Analysis, pp. 295–317. Washington DC: George Mason University Press.
  • Doreian et al. (1984) Doreian, P., K. Teuter, and C.-H. Wang (1984). Network autocorrelation models: some monte carlo results. Sociological Methods & Research 13(2), 155–200.
  • Dow et al. (1982) Dow, M., M. Burton, and D. White (1982). Network autocorrelation: A simulation study of a foundational problem. Social Networks 4, 169–200.
  • Egger and Staub (2016) Egger, P. H. and K. E. Staub (2016). Glm estimation of trade gravity models with fixed effects. Empirical Economics 50(1), 137–175.
  • Fletcher (1970) Fletcher, R. (1970). A new approach to variable metric algorithms. The Computer Journal 13(3), 317–322.
  • Franzese and Hays (2007) Franzese, R. J. and J. C. Hays (2007). Spatial econometric models of cross-sectional interdependence in political science panel and time-series-cross-section data. Political Analysis 15(2), 140–164.
  • Garlaschelli and Loffredo (2005) Garlaschelli, D. and M. I. Loffredo (2005). Structure and evolution of the world trade network. Physica A: Statistical Mechanics and its Applications 355(1), 138–144.
  • Genz et al. (2016) Genz, A., F. Bretz, T. Miwa, X. Mi, F. Leisch, F. Scheipl, and T. Hothorn (2016). mvtnorm: Multivariate normal and t distributions. http://CRAN.R-project.org/package=mvtnorm. R package version 1.0-5.
  • Gleditsch (2013) Gleditsch, K. S. (2013). Distance between capital cities. http://privatewww.essex.ac.uk/~ksg/data-5.html. Accessed: 2017-04-07.
  • Goldfarb (1970) Goldfarb, D. (1970). A family of variable-metric methods derived by variational means. Mathematics of Computation 24(109), 23–26.
  • Haug et al. (2002) Haug, M., M. Langvandslien, L. Lumpe, and N. Marsh (2002, January). Shining a light on small arms exports: The record of state transparency. Occasional Paper 4, Norwegian Initiative of Small Arms Transfers.
  • Hays et al. (2010) Hays, J. C., A. Kachi, and R. J. Franzese (2010). A spatial model incorporating dynamic, endogenous network interdependence: A political science application. Statistical Methodology 7(3), 406–428.
  • Head and Mayer (2014) Head, K. and T. Mayer (2014). Gravity equations: Workhorse, toolkit, and cookbook. In G. Gopinath, E. Helpman, and K. Rogoff (Eds.), Handbook of international economics, Volume 4, pp. 131–195. Amsterdam: Elsevier Science Publishing.
  • Helpman et al. (2008) Helpman, E., M. Melitz, and Y. Rubinstein (2008). Estimating trade flows: Trading partners and trading volumes. The Quarterly Journal of Economics 123(2), 441–487.
  • Henningsen (2013) Henningsen, A. (2013). censReg: Censored regression (tobit) models. https://cran.r-project.org/web/packages/censReg/censReg.pdf. R package version 0.5-26.
  • Holtom et al. (2012) Holtom, P., M. Bromley, and V. Simmel (2012). Measuring international arms transfers. Stockholm International Peace Research Institute.
  • Holtom and Pavesi (2018) Holtom, P. and I. Pavesi (2018). Small arms survey - trade update 2018. http://www.smallarmssurvey.org/fileadmin/docs/S-Trade-Update/SAS-Trade-Update-2018.pdf. Accessed: 2019-06-02.
  • Kang et al. (2013) Kang, L., R. Carter, K. Darcy, J. Kauderer, and S.-Y. Liao (2013). A fast monte carlo em algorithm for estimation in latent class model analysis with an application to assess diagnostic accuracy for cervical neoplasia in women with agc. Journal of applied Statistics 40(12), 2699.
  • Kolaczyk (2009) Kolaczyk, E. D. (2009). Statistical analysis of network data. Methods and models. New York: Springer Science & Business Media.
  • Lacombe (2004) Lacombe, D. J. (2004). Does econometric methodology matter? an analysis of public policy using spatial econometric techniques. Geographical Analysis 36(2), 105–118.
  • Leenders (2002) Leenders, R. T. A. (2002). Modeling social influence through network autocorrelation: constructing the weight matrix. Social Networks 24(1), 21–47.
  • Leifeld et al. (2017) Leifeld, P., S. J. Cranmer, and B. A. Desmarais (2017). tnam: Temporal network autocorrelation models. https://cran.r-project.org/package=tnam. R package version 1.6.5.
  • LeSage and Pace (2008) LeSage, J. P. and R. K. Pace (2008). Spatial econometric modeling of origin-destination flows. Journal of Regional Science 48(5), 941–967.
  • LeSage and Pace (2009) LeSage, J. P. and R. K. Pace (2009). Introduction to spatial econometrics. Boca Raton: CRC Press.
  • Louis (1982) Louis, T. A. (1982). Finding the observed information matrix when using the em algorithm. Journal of the Royal Statistical Society. Series B 44(2), 226–233.
  • Magnus and Neudecker (1988) Magnus, J. R. and H. Neudecker (1988). Matrix differential calculus with applications in statistics and econometrics. Chichester: John Wiley & Sons.
  • Marsh (2017) Marsh, N. (2017). Norwegian initiative on small arms transfers, firearms and ammunition trade data 1992-2014. http://www.nisat.prio.org. Accessed: 27.03.2017.
  • Marsh and McDougal (2016) Marsh, N. J. and T. L. a. McDougal (2016). Illicit small arms prices: Introducing two new datasets. Technical report, Small Arms Data Observatory.
  • Marshall (2017) Marshall, M. G. (2017). Polity IV project: Political regime characteristics and transitions, 1800-2016. http://www.systemicpeace.org/inscrdata.html. Accessed: 2017-06-02.
  • Martinez-Zarzoso and Johannsen (2017) Martinez-Zarzoso, I. and F. Johannsen (2017). The gravity of arms. Defence and Peace Economics, https://doi.org/10.1080/10242694.2017.1324722.
  • Mathai and Provost (1992) Mathai, A. M. and S. B. Provost (1992).

    Quadratic forms in random variables: theory and applications

    London: Taylor & Francis.
  • McLachlan and Krishnan (2007) McLachlan, G. and T. Krishnan (2007). The EM algorithm and extensions. Hoboken: Wiley & Sons.
  • Metz and Ingold (2017) Metz, F. and K. Ingold (2017). Politics of the precautionary principle: assessing actors’ preferences in water protection policy. Policy Sciences 50(4), 721–743.
  • Oakes (1999) Oakes, D. (1999). Direct calculation of the information matrix via the em. Journal of the Royal Statistical Society. Series B 61(2), 479–482.
  • R Core Team (2016) R Core Team (2016). R: A Language and Environment for Statistical Computing. Vienna, Austria: R Foundation for Statistical Computing.
  • Robert and Casella (2004) Robert, C. and G. Casella (2004). Monte Carlo methods. New York: Springer Science & Business Media.
  • Schumacher et al. (2017) Schumacher, F. L., V. H. Lachos, and D. K. Dey (2017). Censored regression models with autoregressive errors: A likelihood-based perspective. Canadian Journal of Statistics 45(4), 375–392.
  • Shanno (1970) Shanno, D. F. (1970). Conditioning of quasi-newton methods for function minimization. Mathematics of Computation 24(111), 647–656.
  • SIPRI (2017a) SIPRI (2017a). Arms transfers database. Accessed: 2017-03-10.
  • SIPRI (2017b) SIPRI (2017b). Arms transfers database - methodology. Accessed: 2017-03-10.
  • Suesse and Zammit-Mangion (2017) Suesse, T. and A. Zammit-Mangion (2017). Computational aspects of the em algorithm for spatial econometric models with missing data. Journal of Statistical Computation and Simulation 87(9), 1767–1786.
  • Tallis (1961) Tallis, G. M. (1961). The moment generating function of the truncated multi-normal distribution. Journal of the Royal Statistical Society. Series B 23(1), 223–229.
  • Thurner et al. (2018) Thurner, P. W., Schmid, C. Christian, Skyler, and G. Kauermann (2018). The network of major conventional weapons transfers 1950-2013. http://journals.sagepub.com/doi/10.1177/0022002718801965. Online First: Journal of Conflict Resolution.
  • Tinbergen (1962) Tinbergen, J. (1962). Shaping the world economy: An analysis of world trade flows. New York Twentieth Century Fund 5(1), 27–30.
  • UCDP (2019) UCDP (2019). UCDP. http://ucdp.uu.se/downloads/. Accessed: 2018-12-01.
  • Vaida and Liu (2009) Vaida, F. and L. Liu (2009). Fast implementation for normal mixed effects models with censored response. Journal of Computational and Graphical Statistics 18(4), 797–817.
  • Ward et al. (2013) Ward, M. D., J. S. Ahlquist, and A. Rozenas (2013). Gravity’s rainbow: A dynamic latent space model for the world trade network. Network Science 1(1), 95–118.
  • Wei and Tanner (1990) Wei, G. C. and M. A. Tanner (1990). A monte carlo implementation of the em algorithm and the poor man’s data augmentation algorithms. Journal of the American Statistical Association 85(411), 699–704.
  • Wilhelm et al. (2012) Wilhelm, S. et al. (2012). Moments calculation for the doubly truncated multivariate normal density. arXiv preprint arXiv:1206.5387.
  • World Bank (2017) World Bank (2017). World bank open data, real GDP. http://data.worldbank.org/. Accessed: 2017-04-01.
  • Youden (1950) Youden, W. J. (1950). Index for rating diagnostic tests. Cancer 3(1), 32–35.
  • Zitzewitz (2012) Zitzewitz, E. (2012). Forensic economics. Journal of Economic Literature 50(3), 731–69.