Estimating Nonlinear Network Data Models with Fixed Effects

This paper considers estimation of a directed network model in which outcomes are driven by dyad-specific variables (such as measures of homophily) as well as unobserved agent-specific parameters that capture degree heterogeneity. I develop a jackknife bias correction to deal with the incidental parameters problem that arises from fixed effect estimation of the model. In contrast to previous proposals, the jackknife approach is easily adaptable to different models and allows for non-binary outcome variables. Additionally, since the jackknife estimates all parameters in the model, including fixed effects, it allows researchers to construct estimates of average effects and counterfactual outcomes. I also show how the jackknife can be used to bias-correct fixed effect averages over functions that depend on multiple nodes, e.g. triads or tetrads in the network. As an example, I implement specification tests for dependence across dyads, such as reciprocity or transitivity. Finally, I demonstrate the usefulness of the estimator in an application to a gravity model for import/export relationships across countries.


Introduction
Networks are common in both economic and social contexts, and it is important to understand the factors that play a role in both the formation and strength of the links between agents. The econometric analysis of networks faces a number of challenges that have received much attention in recent literature (see de Paula (2020) and Graham (2020) for reviews of this literature). One common modeling approach is to assume a dyadic network structure (one in which decisions are made bilaterally between agents), but allow for linking decisions to depend on unobserved agent-specific heterogeneity. These models are common in practice since they are straightforward to implement while still being able to capture important aspects of observed networks. Controlling for agent-specific heterogeneity is important since in many real world networks agents vary significantly in the number and strength of connections made. Ignoring this heterogeneity can lead to large biases in estimated effects.
In this paper, we consider the estimation of dyadic models, where the presence of unobserved heterogeneity is accounted for by two sets of agent-specific fixed effects -a sender and a receiver effect. The fixed effects approach is appealing as it does not require strong assumptions about the unobserved component as in random effects models. In addition, the network setting does not suffer from the 'fixed T ' issues of panel data, since we observe every agent interacting with N − 1 other agents in the network, so that fixed effect estimates are consistent. However, the large number of fixed effects (proportional to the square root of the sample size) does create an incidental parameter problem (Neyman and Scott, 1948). This paper proposes a jackknife approach to bias correction, which has a number of benefits over existing methods. Importantly, the jackknife is easily adaptable to a range of settings, including models for non-binary outcome variables. Additionally, since the jackknife estimates all of the parameters in the model, including the fixed effects, we are able to construct estimates of average effects and counterfactual outcomes. We also show how the jackknife can be used to bias correct averages over functions of multiple observations (e.g. dyads or triads in the network), which we show is useful for constructing specification test statistics, such as tests for the presence of certain strategic interactions like reciprocity or transitivity.
We demonstrate the consistency and asymptotic normality of the jackknife estimator under asymptotic sequences in which a single network grows in size, while the network remains 'dense'. The network model we consider is one in which agents make bilateral decisions about link-specific outcomes, independently of other relationships. This type of dyadic model (a dyad is a pair of agents) has received much attention in the literature, because of its tractability and its ability to replicate some key features of observed networks. In particular, it allows for: homophily, the tendency of agents to form stronger ties with other agents that are similar to them; and degree heterogeneity, where the number/strength of links in the networks can vary substantially across nodes.
In the case of a binary outcome variable, the model we consider is one of link formation, and is an extension of the model by Holland and Leinhardt (1981). There are several alternative approaches to address the incidental parameters problem in this setting. Graham (2017), Charbonneau (2017), and Jochmans (2018) all consider versions of this model in which the latent disturbances follow a logistic distribution, and use conditioning arguments to remove dependence on the fixed effects. The conditioning approach has the advantage of being applicable under certain sparse network asymptotic sequences, but is limited to models in which sufficient statistics for the fixed effects exist, and is not able to recover counterfactuals or average effects. Yan et al. (2019) also studies the logistic model and provides asymptotic results for the incidental parameters. Graham (2017) considers an analytical correction for the logistic model, while Dzemski (2019) derives the analytical correction for a probit model. The analytical bias correction approach is limited to dense network sequences, as in this paper, and similarly to this paper can recover average effects. The advantage of the jackknife correction relative to an analytical approach is that it provides an off-theshelf approach that researchers may apply to new settings, without the need to first derive bias expressions. Candelaria (2020), (Toth, 2017), and Gao (2020) study identification of the common parameters without a known parametric form for the disturbance term, while Zeleneev (2020) allows for nonparametric structure in the unobserved heterogeneity term.
Although the focus of the literature on dyadic network models has been on the binary outcome case, researchers often have access to outcome variables that are non-binary. Examples of these settings include the value of exports between countries, the value of loans between banks, or the number of workers migrating between states. The results in this paper are derived for a general M-estimator satisfying basic regularity conditions and so cover a range of models for both binary and non-binary outcome variables, as well as a range of estimation approaches, including MLE, quasi-MLE and nonlinear least squares estimation.
As a demonstration of the technique in an empirical setting, we estimate a model of international trade relationships. Gravity models have been a workhorse model in the trade literature for many years, and the importance of including country-specific fixed effects is well known ((Anderson and Van Wincoop, 2003)). We estimate the zero-inflated negative binomial model of Burger et al. (2009), which combines both a model for the decision of countries to engage in trade, as well as a model for the value of exports conditional on some trade occurring. The model addresses to key issues in the gravity model literature: it allows for a large proportion of zero trade flows in the network, and it captures the observed high dispersion of export values across countries. We obtain bias-corrected estimates of both the model parameters, as well as average effects.
The jackknife bias correction also allows for the construction of various specification tests. Many models of network formation include strategic aspects in which agents' decisions are influenced by the state of the network. For instance, agent i may find it more beneficial to link with j if they already share many other links in common. Graham and Pelican (2020) derives the locally best similar test for a class of alternatives in a logit model, using conditioning arguments. Dzemski (2019) tests for the presence of transitive links with triads (groups of three agents) in a probit model, and derives an analytical bias correction for the statistic. We demonstrate that a range of test statistics, including that of Dzemski (2019), can be bias-corrected using the jackknife. This extends the set of tests available to researchers, as well as the range of models they can be applied to. As an example, we test for reciprocity and transitivity in trade links between countries and find evidence that the decisions of countries to engage in trade are reciprocal (if country i exports to country j then it is likely that j also exports to i), but do not find evidence of transitive relationships.
The network jackknife extends previous results on jackknife bias correction in panel data. Hahn and Newey (2004) introduced a jackknife correction for panel estimators with individual fixed effects, based on re-estimating the parameters on data sets that exclude a single time period. Dhaene and Jochmans (2015) present a split-sample version of this idea based on estimating the model in the first and second halves of time periods separately. Fernández-Val and Weidner (2016) develop a general framework that allows for both time and individual fixed effects. The analysis in this paper builds heavily off of the asymptotic expansions in Fernández-Val and Weidner (2016).
Analogously to the panel data setting, the network jackknife is constructed by forming 'leaveout' estimates that exclude certain subsets of the data. Cruz-Gonzalez et al. (2017) and Chen et al. (2021) have suggested jackknife approaches for network data, although without formal proof, based on either a split-sample approach or a leave-one-out approach that drops a single agent at a time. We propose a different approach to jackknifing network data that is based on a novel partitioning of the data set that constructs leave-out estimates that remove bias from both sender and receiver of fixed effects in one step. We extend the asymptotic expansions of Fernández-Val and Weidner (2016) to allow for formal analysis of the jackknife estimator. The jackknife proposed here drops a single observation per fixed effect at each step, so that our approach is likely to have better finite sample variance properties than the split-sample approach (see Hughes and Hahn (2020) for a formal argument in the panel setting). In contrast to a jackknife that drops all observations from a single agent in the network, our jackknife retains all agents in each leave-out estimation, so that the distribution of unobserved effects is held constant. We demonstrate the small-sample effectiveness of our approach in comparison to previous suggestions in simulations that show that our jackknife is more robust to settings with meaningful levels of unobserved heterogeneity and in networks with lower density.
In addition, we introduce a weighted jackknife, that differs from standard implementations of the jackknife approach by taking a weighted-average of the leave-out estimates. This version puts less weight on noisier leave-out estimates, which improves the finite-sample properties of the jackknife in sparser settings. The weighted jackknife idea may be useful elsewhere, for example in binary-outcome panel data models with few successes for some individuals (so called 'rare events'). Finally, we also introduce a 'leave-l-out' version of the jackknife. This version requires only (N −1)/l additional estimations of the model, and may allow researchers to reduce the computational burden in settings where model estimation is difficult.
The rest of the paper is organized as follows. Section 2 introduces the network model and discusses implementation of the jackknife procedure for the estimation of model parameters, while Section 3 discusses estimation of average effects, and the construction of specification tests. Section 4 provides asymptotic results for the estimators, and discusses the main assumptions under which they hold. In Section 5 we demonstrate the method by estimating a model of international trade flows, while Section 6 reports simulation results that are consistent with the jackknife theory.

Model
The researcher observes a network of N agents; these agents may, for example, be individuals, firms, or countries. For each potential directed connection, i → j, we observe an associated link-specific outcome variable Y ij . The variable Y ij may capture the presence (or absence) of a link between two agents, in which case Y ij is binary, or may represent a measure of the strength of the link between agents. For example, Y ij may be the value of exports from country i to country j in a particular year, or the number of times agents i and j interacted in some period. Links are directed, meaning that Y ij = Y ji in general, and so, following the literature, we term i the 'sender' and j the 'receiver' in link Y ij .
The researcher also observes a set of link-specific covariates X ij . The covariates capture characteristics of the relationship between agents that may impact the linking outcome. Often these will be interpreted as measures of homophily, that is, the tendency for agents to link with other agents that are similar to themselves. For example, countries may engage in greater levels of trade if they share a common language, or are geographically close.
Agents are endowed with two fixed effects, α i and γ i , which capture unobserved degree heterogeneity, that is, the tendency of some agents to form more (or stronger) links than others. The 'sender' fixed effect α i accounts for heterogeneity in out-degree (the number or strength of links from agent i to other agents), while the 'receiver' fixed effect accounts for in-degree heterogeneity. Degree heterogeneity is an important feature of many networks, for example, we would expect countries with larger GDPs to engage in more trade than smaller countries ceteris paribus (see Anderson and Van Wincoop (2003) for an example of such a model). Since the network considered here is a directed one, we allow for the sender and receiver fixed effects to differ; some countries may have structural tastes for importing goods over exporting, that is, they run trade deficits (or vice versa), so that α i < γ i . Failure to account for degree heterogeneity in a network can lead to incorrect conclusions about the strength of homophily in a network. For example, observing that the United States imports more from China than from Canada may lead to the conclusion that distance between countries is unimportant for trade if we do not account for a China export effect. Graham (2017) provides some further intuition for why failing to account for degree heterogeneity can bias conclusions about homophily in a network.
We make the assumption that linking decisions are bilateral in nature, so that where ⊥ ⊥ denotes independence of the outcomes conditional on observed covariates and fixed effects. This assumption does allow for dependence between the two links within a pair of agents (a dyad), but implies the decision between i and j is independent of that between i and k for instance. Importantly, this independence is conditional on the covariates and agent-specific fixed effects. Unconditionally, country i's exports to country j are correlated with their exports to country k, since both are determined by the exporter effect α i . In many settings, the inclusion of fixed effects will be important in establishing the plausibility of (1). Assumption (1) may not be appropriate in situations where linking decisions are strategic. Estimation of models with strategic interactions is substantially more challenging, and is likely to require multiple observations of the network over time. Nonetheless, the dyadic model presented here still represents an important baseline model, and can be used to construct tests for the presence of strategic interactions against the null hypothesis of (1).
We discuss examples of such tests in Section 3.2.
We leave the specific model for the network outcomes unspecified, and assume only that the model parameters (β 0 , α 0 , γ 0 ) are solutions to the population maximization problem whereĒ represents expectation conditional on the exogenous covariates and fixed effects, and is a known function that is maximized in expectation at the true parameters. Many models can be estimated by maximizing objective functions of the form in (2), including MLE, quasi-MLE, and nonlinear least squares estimators. The researcher need only specify the objective function L and does not need to specify the distribution of the fixed effects, or how they relate to the covariates.
We assume that the unobserved effects enter in an additively separable manner, i.e. as α i + γ j , identification of the two sets of fixed effects parameters requires a normalization, for which we choose is a penalty term intended to impose this normalization on the fixed effect parameters, where b > 0 as an arbitrary constant. 1 Note that the vectors of fixed effects α and γ are dependent on the network size N , although we leave this 1 In practice the constraint could simply be eliminated by substituting it into the objective. We follow Fernández-Val and Weidner (2016) in choosing this normalization as it is convenient in the proofs. dependence implicit in the notation. When functions are evaluated at the true parameters (β 0 , α 0 , γ 0 ) we will typically drop them from notation. We will also use the shorthand ij = (Y ij , X ij , β 0 , α 0,i + γ 0,j ) when convenient.

Example. Maximum likelihood
If we specify the full conditional distribution of the outcome variable as then will be the log-likelihood function Note that L need not be a true log-likelihood, since it may be that the observations Y ij and Y ji are conditionally dependent, in which case L is a pseudo log-likelihood.
As an example, consider a model of directed link formation according to where ε ij follows a known distribution F . This linking rule is compatible with a model in which utility is transferable across linked agents as in Bloch and Jackson (2007). Given the distributional assumption for ε ij , the probability of a link forming, conditional on covariates and fixed effects is p ij = F (X ij β + α i + γ j ), and the log-likelihood is ij = Y ij log p ij + (1 − Y ij ) log(1 − p ij ). This is an extension of the linking model of Holland and Leinhardt (1981) and has been used extensively in empirical literature.

Example. Nonlinear least squares
The researcher may choose to specify only the conditional mean function for the outcome, rather than its full distribution, e.g E Y ij |X, β, α, γ = h(X ij , β, α i + γ j ). In this case, we may estimate the parameters of the model by setting

Empirical setting -gravity equation
Although much of the focus of the incidental parameters bias-correction literature has been on the binary outcome case, researchers often have access to non-binary outcome variables and will be interested in modeling networks in which the links are weighted. As a working example throughout the paper, we will consider a model of country-level trade relationships using a data set consisting of a directed network of export volumes between 136 countries (136 × 135 country pair observations) in 1990. The data are taken from Santos Silva and Tenreyro (2006), and additional details on their construction can be found in their paper. The outcome variable is the value of exports from country i to country j. We also use several covariates to capture homophily in trade relationships, which include: log distance, the log of the distance between the capitals of the countries; border, an indicator of whether the countries share a common border; language, an indicator for whether the countries share a language; colonial, and indicator for whether either country had colonized the other at some point in history; and trade agreement, an indicator for the presence of a joint preferential trade agreement between the two countries.
The Anderson and Van Wincoop (2003) gravity equation expresses the trade volume from country i to country j as where Y ij is the trade flow from country i to country j, G i is GDP of country i, and d ij is inversely proportional to the distance between the two countries (which is generally taken to include all factors that create resistance to trade). The inclusion of exporter and importer fixed effects (α i and γ i ) is intended to control for multilateral resistance terms, which may bias results if excluded.
A simple method for estimating the parameters in (3) is to first log-linearize the model. Unfortunately, this raises the issue of how to deal with the presence of zero outcomes that are common in trade data. In the country-level trade data introduced above, just under half of all country pairs engage in no trade. A number of solutions to this problem have been suggested. Several authors use Tobit models or two-step Heckman style models, which combine a binary selection equation (predicting whether or not any trade occurs) and a separate equation for the value of trade (conditional on selection); see for example Helpman et al. (2008), Rose (2004), Linders and de Groot (2006).
Another popular approach, suggested by Santos Silva and Tenreyro (2006), is to use a Poisson pseudo-maximum-likelihood estimator, which provides a natural way to incorporate zerovalued outcomes, as well as being robust to heteroskedasticity issues that can arise when log-linearizing multiplicative equations. However, there are two key drawbacks to modeling trade flows with a Poisson distribution. Firstly, the proportion of zeroes observed in typical trade data is much larger than that predicted by a Poisson model. Secondly, since the variance of a Poisson is restricted to be equal to its mean, outcomes are typically much more dispersed than would be expected under the Poisson model. In order to address these two issues, Burger et al. (2009) propose a zero-inflated negative binomial model, in which the value of trade between i and j is given by the product of two variables, 1} is a binary decision to enter into a trading relationship, while Y * ij is the value of exports that will be realized, conditional on z ij = 1. The binary decision is modeled using as a probit function, while the latent outcome Y * ij is modeled as a negative binomial variable, which allows for overdispersion in the model for Y * ij , that is, it allows the variance to differ from the mean. Since the distribution of Y ij is parametrically specified, we may estimate the model using maximum likelihood, so that this model is captured by the framework in (2).

Incidental parameters problem
In total, the model contains dim(β)+2N parameters to be estimated, from the N (N −1) observations (Y ij , X ij ) -we will typically refer to these observations using the shorthand (i, j). As is well known, nonlinear estimators with fixed effects suffer from an incidental parameter problem (Neyman and Scott, 1948). To describe the problem, consider the maximization problem (2) after first concentrating out the fixed effect parameters Replacing the population functions α(β), γ(β) = arg max α,γĒ L N (β, α, γ) with their sample values, results in an objective function that is biased, in the sense that To see why, observe that the first-order condition for α i (β) depends only on the N − 1 observations (i, j) for j = i. Similarly, the first-order condition for γ i (β) depends on the N − 1 observations (j, i) for j = i. Expanding α i (β) around α i (β) (and similarly for γ) will therefore result in a bias of order O(N −1 ). Under regularity conditions discussed in Section 4, we show that the bias of the maximizer of the profile objective function (5) is approximately given by for some bias term B N . Analogously to the panel literature, the bias is inversely proportional to the number of observations used to estimate each of the fixed effect parameters; the exporter effect for country i is estimated using the data on the N − 1 other countries in the network that i may export to. As the size of the network grows, N → ∞, we will have that β N → β 0 so that parameter estimates are consistent. Considering β as an estimator for β N , we can show that N ( β − β N ) ⇒ N (0, V ). However, since the bias β N − β 0 is of the same order as the estimation error, O(N −1 ), β will be asymptotically biased, that is The incidental parameters generate an asymptotic bias in the network model analogous to the panel setting with both N and T growing to infinity at the same rate. Similar asymptotic expansion arguments have been used in the panel data literature on nonlinear models with fixed effects. Hahn and Newey (2004) derive expansions for models with individual fixed effects, while Fernández-Val and Weidner (2016) derive expansions that apply to general models with additively separable unobserved effects, and Chen et al. (2021) consider the setting with interactive effects. The expansions used in this paper rely heavily on these prior results. Dzemski (2019) also applies the Fernández-Val and Weidner (2016) expansions to the network model structure to derive bias expressions for a probit model. In this paper, we extend the asymptotic expansions to higher order so that they may be used to justify jackknife bias correction procedures. We demonstrate that, under mild additional regularity conditions, the jackknife estimator introduced below is asymptotically normal and mean zero, so that valid inference can be performed on model parameters.

Jackknife bias correction
The jackknife bias-corrected estimator is constructed as a linear combination of the fullsample parameter estimates and an average of 'leave-out' estimators that exclude certain observations in the data set. The particular linear combination chosen can be motivated by asymptotic expansions of the estimator, in particular the form of the bias in (6).
Suppose we were to drop observations from our data set in such a way that for every country i we exclude one observation in which i is the exporter, and one observation in which i is the importer (recall that a single observation is one export-import relationship, of which we observe N (N − 1) in total). We show that the new estimator using this 'leave-out' sample, β, has a bias that is approximatelyĒ The form of the bias forβ can be explained by two important factors: (i)β is estimated using only N − 2 observations per fixed effect, since we excluded one observation related to each fixed effect parameter, and (ii)β is estimated on a random sample generated from the same set of fixed effects, so that the bias expression B N is the same as that in (6).
Taking advantage of the fact that the estimateβ has a larger bias than the full-sample estimate β by the factor N −1 N −2 , we can construct a new estimator β jack = (N − 1) β − (N − 2)β which has no asymptotic biasĒ To describe the construction of the leave-out estimators, we first define a partition of the N (N − 1) observations of directed pairs (i, j) into N − 1 sets of the form 2 I k = (i, j) : j = (i + k) (mod N ) , for k = 1, . . . , N − 1 that is, the set of directed pairs {(1, 1 + k), . . . , (N − k, N ), (N − k + 1, 1), . . . , (N, k)}. 3 Figure 1 represents the structure of the first three sets I k for a network of N = 10 agents. Observations are ordered in an N × N matrix so that the (i, j) cell represents the corresponding observation in the network (the diagonal elements are empty since there are no (i, i) observations). The leave-out sets take diagonal sections from the data matrix. Importantly, constructing the sets this way ensures that each contains exactly one observation related to each sender and receiver fixed effect; i.e. there is one observation taken from every row and every column. Let 1 k ij = 1{(i, j) ∈ I k }, be an indicator variable that is equal to one whenever the observation (i, j) is not in the k-th leave-out set. The k-th leave-out estimates are subject to that is, the estimates obtained by excluding the observations in I k from the data. We can then construct the jackknife bias-corrected estimator The construction of the leave-out estimators is analogous to jackknife bias correction in the panel data setting; however, the structure of the jackknife proposed here is new. The procedure relies on dropping sets of observations that contain a single observation related to every sender fixed effect α i as well as every receiver fixed effect γ i . In this way, the bias from both types of fixed effects can be addressed simultaneously, while holding the distribution of fixed effects constant across the leave-out samples. This is in contrast to an approach which drops all observations from a single agent, which removes that agents' fixed effects from the leave-out sample and alters the distribution of unobserved heterogeneity. We show in simulations that the method proposed here is more robust to networks that have more unobserved heterogeneity or are less dense.
We prove in Section 4 that the jackknife bias correction is consistent, and asymptotically normal, with mean zero and variance equal to that of the full-sample estimate The fact that the jackknife is able to remove bias without affecting the asymptotic variance of the estimator may seem surprising. This important feature is achieved by averaging across the N − 1 different leave-out estimators β (k) . Since the sets I 1 , . . . , I N −1 form a partition of the N (N − 1) observations in the network, each observation is excluded from exactly one of the leave-out estimates. This balanced treatment of observations ensures that the jackknife procedure does not affect the first-order asymptotic approximation of the estimator. 4 Remark 1. The construction of the leave-out sets depends on the labelling of nodes, and so the final estimator will be dependent on this labelling (since the make-up of the leave-out sets will change). While the researcher could re-randomize the node labels, construct β J for each randomization, and then average to remove some of the arbitrariness of the node labels, this is not necessary. The estimators with different labelings should be very similar so that the additional computations will have little effect. 5 The jackknife estimator requires N estimations of the model, and so may be computationally intensive for large networks, although speed may be improved by computing the leave-out estimates in parallel, and using good starting values such as the full sample estimates. As an alternative, we present a 'leave-l-out' version of the jackknife, which reduces the number of additional estimations of the model by dropping l observations per fixed effect, as opposed to just one. To describe the estimator, let N l = N −1 l (we assume here that N −1 is divisible by l for simplicity). We can construct the k-th leave-l-out set by combining l of the leave-one-out sets as follows I l k = ∪ l−1 j=0 I k+jN l , for k = 1, . . . , N l . This results in N l = N −1 l non-overlapping leave-l-out sets, with corresponding estimatesβ l,(k) , which are the estimates from using all 4 In the panel data setting, Dhaene and Jochmans (2015) note that forming a jackknife using overlapping subpanels (across time) results in an inflation of the asymptotic variance, since some time periods are used more that others. That the I k form a partition ensures that each observation is used an equal number of times in the average 1 N −1 N −1 k=1 β (k) . 5 In simulations (see Section 6 for the design) the estimates based on different labelings were close to identical and so the additional calculations had almost no effect on the estimation. observations except those in the k-th leave-l-out set I l k . A jackknife bias-corrected estimate can then be constructed aŝ where N l = N −1 l . 6 Remark 2. The leave-l-out jackknife bias correction has the same asymptotic variance as the standard leave-one-out jackknife and the full-sample estimator. However, there may be some finite-sample efficiency loss, particularly when l is large or when the network is not sufficiently dense. Hughes and Hahn (2020) show in the panel case that the leave-one-out jackknife is higher-order more efficient than the split-sample jackknife (i.e. its variance to O(N −1 ) is smaller), and it is likely that the same result applies here, although this is beyond the scope of the present paper.

A weighted jackknife
The jackknife relies on large dense network asymptotics, but in finite samples it is possible for some leave-out estimates to drop a number of important observations all at once. This is more likely to occur when N is small, there are few links for some nodes, or when we are using the leave-l-out jackknife with large l.
The performance of the jackknife can be improved in these settings by taking a weighted average of the estimates β (k) . Define the weights . 6 In the case that N − 1 is not divisible by l, let N l = N −1 l . Then we may partition the data into r = N − 1 − lN l leave-(l + 1)-out sets and N l − r leave-l-out sets. Denoteβ (k) as the leave-(l + 1)-out estimates for k = 1, . . . , r and the leave-l-out estimates for k = r + 1, . . . , N l . Then, let The jackknife bias-corrected estimator is then given bŷ The weighted-jackknife estimator is given by . The weights W (k) are the Hessian for β, after concentrating out the fixed effects. In the special case that L is a log-likelihood function, W N is the Fisher information for β, and so is equal to the inverse of the asymptotic variance. In this case, we are using an inverse variance weighting scheme, which down-weights leave-out samples that produce particularly noisy estimates of the common parameters. The weighting scheme is equally applicable to non-likelihood settings, although it no longer carries the inverse variance interpretation. In simulations, this weighted version of the jackknife significantly improves the performance of the estimator in sparser networks (see Section 6 for more details).
Remark 3. Asymptotically the weights have no effect on the estimator, since all W (k) converge to the same quantity. This implies that the asymptotic variance of the weighted jackknife is the same as that of the standard jackknife. In finite samples, variation in the weights depends on the number of nodes N , as well as the density of the network (i.e. the variation in outcomes for each node). The weighting scheme is likely to have a large impact for small or less dense networks, but in denser (or larger) networks we will havē W −1 J W (k) ≈ I dim β , so that the weighted and unweighted jackknife estimates are very similar.
Remark 4. The motivation behind the particular choice of weights comes from the firstorder asymptotic expansion of the estimator. A Taylor expansion of the first-order conditions of the objective function gives an expression of the form where A is mean zero and asymptotically normal with varianceΩ (as in Theorem 1), while B is an additional term responsible for the asymptotic bias of the estimator.
As demonstrated in Lemmas 3 and 4 in the Appendix, the jackknife procedure applied to A + B successfully demeans B and leaves A unchanged, i.e.
This results in an estimator with no asymptotic bias and unchanged asymptotic variance. In practice however, the jackknife is applied to β, so we must consider the effect of the jackknife procedure on W −1 N (A + B). The validity of the jackknife relies on the fact that W (k) ≈ W N , which is guaranteed in large samples under our assumptions, in particular, dense network asymptotics. However, in finite samples, W (k) could vary substantially in some leave-out samples when N is small, there are few links for some nodes, or we are using the leave-l-out jackknife with large l. This motivates instead averaging over W (k) β (k) , so that variation in W −1 (k) has less impact on the quality of the asymptotic approximation. Intuitively, we are jackknifing the first-order condition for β, rather than β itself. 7

Estimating average effects
In addition to estimation of the common parameter β, researchers may also be interested in estimating certain averages over the distribution of exogenous regressors and fixed effects. An important advantage of the jackknife bias correction, over methods based on conditioning on sufficient statistics (e.g. Graham (2017), Jochmans (2018)), is that by estimating the fixed effect parameters we are able to construct estimators for these averages. Common examples include average and marginal effects, as well as counterfactual outcomes. In the network setting, these are averages over functions of a single potential link (i, j) in the network. We will additionally show that averages over multiple links also provide interesting objects of interest; for example, averages over dyads {(i, j), (j, i)}, triads (groups of three nodes), or other network patterns. As an example, we focus on how these objects can be used to construct tests of the assumption of independent link formation stated in (1), but they may have wider relevance in empirical work. Estimation of the many fixed effect parameters means that these averages also suffer from an incidental parameter problem. We show that the jackknife can be used to bias-correct average effects estimates and obtain correct inference.

Averages over single observations
A simple fixed effect average may be expressed as where the expectation is taken over the joint distribution of covariates X ij and fixed effects (α i , γ j ), and the function m represents the effect of interest. Here we specify two possible parameters of interest, δ the population average, and ∆ N , the sample average; this choice will affect the asymptotic distribution of the estimator, a point we return to in Section 4. As earlier, we will impose that the fixed effects enter the function m in an additively separable way, as π ij = α i + γ j ; this will imply that the choice of fixed effect normalization will not affect the estimator.

Example. Marginal effect
As an example, consider a binary outcome model with P (Y ij = 1|X, β, α, γ) = F (βX ij + α i + γ j ). We may be interested in estimated the average partial effect of the covariate, in which case we would have Alternatively, we may be interested in the average partial effect at some fixed value of

Example. Counterfactual change
Alternatively, assume that we estimate the conditional mean function E[Y ij |X ij , β, φ] = h(βX ij + α i + γ j ). We may be interested in the counterfactual change in predicted outcome from a change in the value of the covariate X ij from x 0 to x 1 , e.g. the effect of entering or exiting a trade agreement. In this case The average effect in (11) can be estimated by plugging in estimates of the model parameters As with estimates of the parameters themselves, the average effect estimate ∆ N is asymptotically biased, that is, N ( ∆ N − ∆ N ) converges to a normal distribution that is not centered at zero. The asymptotic bias in ∆ N stems from three sources: (i) bias in the common parameter estimates β, (ii) averaging over a nonlinear function of noisy fixed effect estimates (a Jensen inequality type bias), and (iii) correlation between the fixed effect errors and m(X ij , β, α i , γ j ).
The average effect estimator can be bias-corrected using the jackknife in an almost identical way to the bias correction of β. Let be the average effect estimate that uses the parameter estimates in the k-th leave-out sample (that is, the sample which we drop observations (i, j) ∈ I k ). A jackknife bias-corrected estimator is We show in Section 4 that the bias-corrected estimator is asymptotically normal with mean zero. Note that there is no need to use a bias-corrected estimate of β in the construction of the average effects. The jackknife takes care of the bias generated by bias in β as well as the other sources of bias in a single step.

Specification testing
The parameter in (11) is an average over a function that depends on a single observation in the network. In some cases, we may be interested in averages over functions that depend on multiple observations such as patterns depending on pairs of nodes (dyads), groups of three or four nodes (triads or tetrads), or other structures. Although these averages may be of interest in their own right, they prove particularly useful in developing specification tests, and we focus on this case. In this section, we show that like simple average effects, these averages also suffer from the incidental parameters bias problem, but can be bias corrected using the jackknife approach.
Let λ be a set of observations in the network; for example, λ = {(i, j), (j, i)} collects the two observations within a dyad, and λ = {(i, j), (j, k), (k, i)} collects a sequence of potential links between three nodes. Let Λ N be the set of all possible λ formed by permuting the nodes for a network of size N . We consider averages of the form where ∈λ collect the outcomes, covariates and fixed effects for the observations in λ. These generalize the averages in (11) in two ways: (i) they allow for averages over functions of multiple observations in the network, and (ii) they allow the function m to depend on the outcome variable Y ij .
One important application of the type of averages in (13) is to specification testing. The model presented in this paper assumes that decisions are made bilaterally, that is, agents i and j decide on Y ij independently of other outcomes in the network. In some settings, we may expect that decision making has some strategic aspect, in that an agent's utility from a link depends on the presence (or strength) of other links. One way to model such a phenomenon is to include network statistics in the utility function. For example, imagine that i sends a link to j according to where S ij is the value of some network statistic, and ε ij are independent shocks. Models of this form generally result in multiple equilibria, which raises a number of difficulties for estimation. However, under the null hypothesis H 0 : δ = 0, the model is the dyadic link formation model considered in this paper and can be consistently estimated. This suggests that a test statistic based on S ij may be useful for testing the null hypothesis of the dyadic model (i.e. assumption (1)) against an alternative of the form (14). One possible test statistic is where The statistic tests the ability of the network statistic S ij to predict the model errors and so values of T N far from zero suggest the presence of unexplained strategic interactions. Graham and Pelican (2020) consider a similar setup in the case of a logit model and derive the locally best conditional similar test for the null hypothesis H 0 : δ = 0 in (14). The resulting statistic is similar to (15), although this certainly does not imply any optimality of the test proposed here. The motivation for the statistic suggested here is heuristic, but has the advantage of being applicable to a wide set of models (e.g. models that do not admit a sufficient statistic) and has the simplicity of using asymptotic critical values. Some examples of potential test statistics in this framework may be useful.

Example. Reciprocity in link formation
Consider a model in which links are reciprocal, that is the presence of a link from j to i increases the utility of the reverse link from i to j. In this case we could let S ij = Y ji in (15), which gives the statistic This statistic measures whether the average prediction error for reciprocal links differs from the average prediction error for non-reciprocal links. Note that reciprocity is allowed for under the assumptions of this paper, so that a rejection of the null hypothesis of no reciprocity does not affect the interpretation of model estimates.

Example. Transitivity
Consider a triad, three nodes in the network (i, j, k). Under the dyadic model, linking decisions are independent across the three pairs of nodes in the triad, but in many settings it may be reasonable to think that the existence of links between two of these pairs may affect the formation of the links in the third. For example, imagine that i has formed a directed link to j, and j has formed a link to k. We may expect the existence of the indirect path from i to k (passing through j) to increase the likelihood of observing the direct link from i to k. This linking structure, i → j, i → k,k → j (shown in the right diagram of Figure 2), is known as a transitive triangle.

Figure 2: Transitive triangles
There are six potential transitive triangles that may exist within any triad; the figure above shows two of these (in which i is the sender of two links). Additional links may exist within the triadthese do not affect the existence of a transitive triangle.
Transitivity in linking is a feature of many models of strategic network formation and we may test the null hypothesis that linking decisions are dyadic against an alternative in which transitivity exists by choosing S ij = 1 N −2 k ={i,j} Y ik Y kj and using the statistic Note that in both examples, and more generally, the outcome Y ij need not be binary. The test applies equally to non-binary outcomes, for example, in a trade network the presence of large export flows of a particular good from i to k and from k to j may reduce the expected direct exports of that good from i to j.
The framework in (14) generates just one possible set of specification tests, and many others statistics of the form (13) are possible. One alternative method is to compare the observed frequency of some possible subgraph configuration with the expected frequency under the assumed dyadic model. Such a test, for the case of transitive and cyclic triangles, was proposed by Dzemski (2019), who also derived an analytical bias correction for the statistic in a binary outcome model with normal disturbances. The statistic suggested by Dzemski (2019) is of the form where Y ij Y ik Y kj is an indicator for a transitive triangle and p ij p ik p kj is the probability of observing such a triangle when all three links are independent. Many test statistics of this form could be derived in the same way, by taking some function of network outcomes between multiple agents and comparing it to its expectation under the dyadic model. For example, an alternative statistic for reciprocity would be A key advantage of the jackknife procedure proposed in this paper is that it applies to such a wide variety of statistics. Given the plethora of potential test statistics, an analysis of their power properties would certainly be useful, although is beyond the scope of this paper.
Like estimates of the common parameter β, the test statistics discussed above suffer from an incidental parameter bias. Although the infeasible test statistics are mean zero under a correctly specified model, the feasible versions, which replace parameters β 0 , α 0 , γ 0 with their estimated values, have an asymptotic bias that leads to incorrect inference. However, we show in Section 4 that statistics of the form in (13), may be jackknife bias corrected.
To describe the jackknife bias correction for these statistics, denote the number of observations contained in λ as r. Let 1 k λ = (i,j)∈λ 1 k ij be an indicator that is zero whenever any of the observations in λ are included in the k-th leave-out set I k . Define the leave-out estimate where β (k) , π λ,(k) are parameter estimates from the k-th leave-out estimation. The factor N −1 N −r−1 accounts for the fact that m λ is dropped from the average whenever any of the r observations it is a function of are dropped. 8 A jackknife bias-corrected estimator is again given by In Section 4, we show that the bias-corrected statistic is asymptotically normal and well centered. In the case of the specification test statistics discussed above, we haveĒ[ where the form of the variance is shown in Theorem 2. This allows us to test hypotheses in the usual way, comparing N ∆ J / √ V ∆ to the quantiles of a standard normal distribution. Further details on the implementation of these tests are discussed below.

Asymptotic theory
We consider an asymptotic framework in which a single network of N agents grows in size, i.e. N → ∞. Recall that the parameters of interest maximize the objective function in (2) for some function ij = (Y ij , X ij , β, α i + γ j ) of the observables Z ij = (Y ij , X ij ) and additive unobserved fixed effects α i + γ j . The asymptotic theory relies on expansions of the objective function, for which we require certain differentiability and moment conditions. We let φ = (α , γ ) denote the 2N × 1 vector of fixed effect parameters and π ij = α i + γ j represent the additive index through which they enter the objective function. We denote derivatives of the function with respect to parameters by ∂ β ij (β, α, γ) = ∂ ij (β, α, γ)/∂β, When evaluating these objects at the true parameter values, we simply write ∂ π q ij and so on.

Asymptotic analysis for the common parameters
The results below are derived under the following set of assumptions. Proofs are provided in the Appendix.
(i) Conditional on (X, α, γ), dyads are independent, that is, For all N , the objective function L is strictly concave over R dim β+2N , and the matrixH = −∂ φφ L is positive definite.
(iii) For all (i, j), the function (β, π) → ij (β, π) is five times continuously differentiable over B ε,ij almost surely. For all (i, j), the partial derivatives of ij with respect to the elements of (β, π) up to fifth order are bounded in absolute value uniformly over (β, π) ∈ B ε,ij by a function M (Z ij ) > 0 a.s., where max i,jĒ [M (Z ij ) 16 ] is a.s. uniformly bounded over N .
The limits plim N →∞W N = W and plim N →∞Ω N = Ω exist for W and Ω positive definite matrices.
(v) There exist constants b min and b max such that for all (β, π) ∈ B ε,ij a.s. uniformly over i, j and N .
Assumption 1(i) specifies that outcomes depend on dyad specific variables only, and not on other features of the network. Conditional on the observed covariates and fixed effects, the outcome Y ij is independent of other outcomes in the network, with the exception of Y ji . Note that unconditionally outcomes are allowed to be dependent, through dependence across covariates X ij and the fixed effects.
Assumption (ii) contains the parametric restriction of the model and requires that the true parameters β 0 , α 0 , γ 0 are solutions to the first-order equations of the objective function. Concavity of the objective function ensures that the population problem has a unique solution. This is satisfied in many common nonlinear models, including the class of regression models with log-concave densities (as well as censored and truncated versions of these models), which includes probit, logit, ordered probit, Tobit, gamma and beta models among others (see Pratt (1981), Newey and McFadden (1994)).
Assumption 1(iii) provides basic smoothness conditions for the objective function. The derivative and moment conditions are required to ensure the validity of the asymptotic expansions to high enough order to establish the properties of the jackknife procedure, and to ensure that remainder terms are well bounded. Analysis of the jackknife requires higher order expansions than are required for characterization of the analytical bias and first-order asymptotic properties of the estimator, and so Assumption 1(iii) is somewhat stronger than the equivalent assumption employed in Fernández-Val and Weidner (2016).
Assumption 1(iv) ensures that the variance of β is non-degenerate. The termW N is the Hessian matrix for the common parameters β, after concentrating out the fixed effect pa-rameters, whileΩ N is the variance of the score for β (again after concentrating out the fixed effect parameters). To describe estimators of these terms, let . This term is the score for β after partialling out the fixed effect parameters. Estimators of the variance terms can be created in the usual way by plugging in estimates of the model parameters, i.e. , where the terms ∂ ββ L, D β ij etc. are evaluated at the estimates β, α, γ. Note that the estimator Ω N allows for correlation between the Y ij and Y ji outcomes.
Finally, Assumption 1(v) ensures that the Hessian for the fixed effect parameters is positive definite. This requires sufficient variation in the outcomes across both dimensions -i.e., variation in Y ij over j (for fixed sender i) and over i (for fixed receiver j). In the binary outcome case, it implies that the model generates a dense network, one in which the number of links formed by each node tends to infinity as the size of the network grows. The assumption may not be reasonable in all empirical settings -in simulations we investigate the robustness of the estimator to sparsity in finite samples. The density assumption can be avoided in settings where sufficient statistics for the incidental parameters exists, such as the conditional logit framework, since estimation of the fixed effects is avoided. This comes at the expense of no longer being able to estimate average effects or counterfactual outcomes.
We now state the main theorem of the paper, on the asymptotic distribution of the jackknife bias-corrected estimator.
Theorem 1. Let Assumption 1 hold and let β J be the jackknife bias-corrected estimator (8), the leave-l-out estimator (9) or the weighted jackknife (10). Let V N =W −1 NΩ NW −1 N and assume that V = plim N →∞ V N exists and is positive definite. Then, N be an estimator of the asymptotic variance, where W N and Ω N are the plug-in estimators shown in (19).
The jackknife estimator is asymptotically normally distributed and unbiased. It also has the same asymptotic variance as the non-bias-corrected estimator in (2). The variance is the usual sandwich form one, and is easily computed. In the case of maximum likelihood we will have thatW N =Ω N so that the variance simplifies to V N =W −1 N . In general this will not be true, for example the researcher may wish to allow for correlation between Y ij and Y ji by clusteringΩ N at the dyad level as in (19).

Asymptotic analysis for fixed effect averages
Here we present asymptotic results for averages of functions that may take more than one observation as arguments. This structure will cover a number of interesting cases such as standard average effects, averages over dyads, triads or other structures in the network, and specification tests. Recall that λ is a set of observations (i, j), and Λ N is the collection of all such sets formed by permuting the nodes in λ. We let Y λ = {Y ij } (i,j)∈λ , X λ = {X ij } (i,j)∈λ , and π λ = {α i + γ j } (i,j)∈λ collect the outcomes, covariates and fixed effects for the observations in λ. The function of interest is m λ = m(Y λ , X λ , β, π λ ), which is a function of (Y ij , X ij , α i + γ j ) for each (i, j) ∈ λ.
It will be useful to define three separate quantities Here, ∆ N is the average effect computed in the observed sample (at the true parameter values),∆ N is the expectation of the average effects conditional on the distribution of covariates and fixed effects in the observed sample, while δ N is the population expectation. We can decompose the estimation error ∆ N − δ into three sources The first term, ∆ N − ∆ N , represents variation caused by estimation of the parameters in the model, including fixed effects. The next term, ∆ N −∆ N , is variation of the sample outcomes m λ around their conditional expectationsm λ =Ē[m λ ]. In the case that m is a function of the data only through X λ , i.e. it does not depend on outcomes Y λ , we will havem λ = m λ and this second term will vanish. Finally,∆ N − δ captures differences in the distribution of covariates and fixed effects in the observed network, relative to the population. In the case that the full network is observed, or whenever∆ N = 0 as is the case for specification tests, we will have that∆ N = δ N .
The results below will rely on Assumption 1, as well as additional restrictions on the choices of λ and m.
Assumption 2. Let λ be a set of r observations (i, j) containing p distinct agents. For every λ, let B ε,λ be a subset of R dim β+r that contains an ε-neighborhood of (β 0 , π 0,λ ) for all N , with ε > 0.
is five times continuously differentiable over B ε,λ a.s. For all λ, the partial derivatives of m with respect to the elements of (β, π λ ) up to fifth order are bounded in absolute value uniformly over (β, π λ ) ∈ B ε,λ by a function M (Z λ ) > 0 a.s., and max λĒ [M (Z λ ) 16 ] is a.s. uniformly bounded over N .
Assumption 2 (i) restricts m to be a function of a fixed number of edges in the network, which allows us to construct leave-out sets to bias correct using the jackknife. It also assumes that ∆ N is an average of all possible arrangements of the nodes in λ, ensuring that averaging occurs over all dimensions. For example, an average of the form 1 N j =i m(X ij , α i , γ j ) is not allowed since we are only averaging over the receiver dimension, while holding i fixed.
Assumption 2 (ii) is analogous to Assumption 1 (iii), and imposes the same differentiability and moment conditions on m as are imposed on . This allows for asymptotic expansions of ∆ N to be derived, in the same way as for β. Finally, (iii) ensures that the unconditional second moments of m are well defined.
The asymptotic distribution depends on the choice of target parameter, either a conditional or population average. The following theorem states the asymptotic result for the jackknife bias-corrected estimator of the conditional fixed effect average∆ N .
Theorem 2. Let Assumptions 1 and 2 hold, and let ∆ J be the jackknife bias-corrected estimator (18). Then If we additionally assume that E (m λ −m λ )(m λ −m λ ) = 0 for sets λ and λ that share exactly one observation in common, then the asymptotic variance is If we have either m λ =m λ or E (m λ −m λ )(m λ −m λ ) = 0 for sets λ and λ that share exactly one observation in common, then the asymptotic variance is In either case, let V ∆ be the plug-in estimator for V ∆ that replaces the unknown θ with estimates θ. Then V ∆ → V ∆ .
Some explanation for the form of the variance may be useful. The two terms h ij and s ij relate to the first two components of (20). The first component ∆ J − ∆ N contains variation from estimation of the common parameters β and fixed effects α, γ. In the Appendix it is shown that this term can be approximated using the delta-method Note that, replacing the jackknife estimate ∆ J with the standard estimator ∆ N would result in additional terms appearing in the above first-order approximation, related to the incidental parameter bias. The second component is The variance of this term depends on the conditional covariance between m λ and m λ for distinct sets λ and λ . Note that if λ and λ share no dyads in common then they are conditionally independent. The variance of this term depends on the condition E (m λ − m λ )(m λ −m λ ) = 0 for sets λ and λ share exactly one observation in common. Under this condition, the variance of λ (m λ −m λ ) is dominated by covariances between m λ and m λ for λ and λ that share exactly one dyad in common; although λ with two or more common dyads also contribute to the variance, there are an order of magnitude fewer such combinations, and so these represent smaller order contributions that do not appear in the asymptotic variance. In settings where E (m λ −m λ )(m λ −m λ ) = 0 for λ and λ that share exactly one observation in common, ∆ N −∆ N is a degenerate U-statistic and its variance is asymptotically of smaller order than the variance from parameter estimation, i.e. ∆ J − ∆ N , and so may be ignored.
Theorem 2 shows how we may construct confidence sets for the parameter of interest∆ N . When the object of interest is the unconditional average δ, the convergence of the estimator will be dominated by variation from the third component in (20),∆ N − δ. To describe the statistic in this setting, it is useful to use its U-statistic representation. We will additionally assume that X ij = h(X i , X j ). This condition appears in other work on dyadic models, for example in Graham (2017). When X ij measures the similarity (or difference) between i and j in some measure, we will commonly have X ij = d(X i − X j ), for d some distance function. Alternatively, if X ij captures common membership in some group, we may have X ij = X i X j where X i is an indicator for i's membership.
To give∆ N a U-statistic representation, we first sum together allm λ which share the same set of agents. Since there are p agents in each λ, this gives p! different λ that can be created from a given set of agents. We denote these sets of unordered agents by τ . We have, for where u τ = 1 p! λ∈τm λ . The term u τ is a symmetric function of {β, X i , α i , γ i } for p agents i. Assuming that the {X i , α i , γ i } are i.i.d. over agents,∆ N − δ is a U-statistic of order p and we may apply standard theory on such statistics to compute its asymptotic distribution.
Theorem 3. Let Assumptions 1 and 2 hold. Additionally, assume that where X i is an observed agent-specific characteristic, and also that ( Additionally, the variance estimator V δ in (21) is consistent, i.e. V δ → V δ .
The convergence rate in Theorem 3 is slower than the rate in Theorem 2. While m λ and m λ are conditionally independent when λ and λ share no dyads in common, the two are unconditionally independent only when they share no agents in common. Since there are many more sets λ that share a single agent i than share a dyad (i, j), the variance of∆ N − δ is an order of magnitude larger than that of ∆ N − ∆ N , and so the convergence rate is slower.
Similarly to Theorem 2, the variance is dominated by covariances between sets λ that share exactly one node in common. The term Σ 1 is the covariance between u τ and u τ when τ and τ share exactly one agent in common. A consistent estimator for this quantity is given by whereμ i is the average over all sets λ containing agent i, m λ is a plug-in estimator form λ , and µ = λ m λ is the overall mean.
Since the rate of convergence in Theorem 3 is N −1/2 , there is in fact no asymptotic bias generated by the incidental parameters. The bias from the estimation of the fixed effects is of order N −1 , which is smaller than the variation in the sampled distribution of fixed effects around is population distribution. Nonetheless, we would still recommend bias correcting estimators as it is likely to improve the finite sample properties of inference, in terms of correct centering of confidence sets, with little or no cost in terms of additional variance.
In the panel data setting, Fernández-Val and Weidner (2016) report such improvements in simulations.

Examples
Theorems 2 and 3 suggest that the construction of confidence sets for estimates and hypothesis testing may be performed in the usual way, by using the asymptotic normal approximations. Standard plug-in estimates for the variance expressions may be used. Here we provide some examples of how these results can be used.

Example 1. Average marginal effect in a probit model
The average marginal effect in the probit model can be estimated using where ϕ is the standard normal density function. In this setting we have λ = (i, j) and m ij =m ij so that s ij = 0. The variance in Theorem 2 is then the standard delta-method variance of and standard plug-ing estimators may be used. A (1 − α)-per cent confidence set could be constructed for ∆ N using ∆ J ± c 1−α/2 V 1/2 ∆ /N , where c 1−α/2 is the (1 − α/2) quantile of the standard normal distribution.
If instead we are interested in inference with respect to the population parameter δ, Theorem 3 states that we may compute the asymptotic variance as

Example 2. Testing transitivity in a probit model
Recall that a statistic for testing the presence of transitivity is Fitting Y ij with a probit regression, we have p ij = Φ(β X ij + α i + γ j ). Sincem λ =Ē (Y ij − p ij )Y ik Y kj = 0, we have∆ N = δ = 0 so we may determine the asymptotic distribution of the test statistic using Theorem 2. We have

From the likelihood for a probit model, we have
Also, from which we can construct h ij using the estimated Hessian matrix and the formula given in Theorem 2.

Empirical example
We illustrate the jackknife procedure on a data set consisting of a directed network of export volumes between 136 countries (136 × 135 country pair observations) in 1990. The data are taken from Santos Silva and Tenreyro (2006), and additional details on their construction can be found in their paper. The outcome variable Y ij is the value of exports from country i to country j. We also use several covariates to capture homophily in trade relationships, which include: log distance, the log of the distance between the capitals of the countries; border, an indicator of whether the countries share a common border; language, an indicator for whether the countries share a language; colonial, and indicator for whether either country had colonized the other at some point in history; and trade agreement, an indicator for the presence of a joint preferential trade agreement between the two countries. Burger et al. (2009) propose a zero-inflated negative binomial model. The value of trade between i and j is given by the product of two variables Y ij = z ij Y * ij , where z ij ∈ {0, 1} is a binary decision to enter into a trading relationship, while Y * ij is the value of exports that will be realized, conditional on z ij = 1. The binary decision is modeled using as a probit function, while the latent outcome Y * ij is modeled as a negative binomial variable.

Zero-inflated binomial model
In this example, the objective function is given by where θ = (β, α, γ, ν), and The parameter ν captures the degree of overdispersion in the model for Y * ij , with ν → ∞ resulting in a model with equal mean and variance (as in the Poisson), while smaller ν lead to greater degrees of dispersion.
Estimates of the parameters in the model are presented in Table 1. Most variables change by only small amounts after bias correction. However, the effect of sharing a common border on the probability of engaging in zero trade changes significantly after bias correction; while the maximum likelihood estimate suggests that common borders are important for link formation, the bias corrected estimate is no longer significant. This suggests that the sharing a common border has little effect on the likelihood of engaging in trade, but does affect the volume of trade. The results also suggests a substantial impact of trade agreements, both on the probability of engaging in trade and on the volume of trade, a result that is robust to bias correction. The overdispersion parameter ν is less than a half, suggesting a significant amount of overdispersion, i.e. export volumes have far greater variation across country pairs than suggested by a Poisson model.
The rightmost column in the table reports the difference between the MLE and jackknife bias-corrected estimators in terms of their standard errors. For a number of variables in the model of export volumes, the change in estimate is around three-quarters of the standard deviation or more, which has an important impact on inference. To give some idea of the scale of these biases, a bias of three-quarters of a standard error is enough for a five per cent test two reject around 12 per cent of the time (more than twice nominal size), while bias of 1.5 standard errors leads to a rejection rate of more than 30 per cent.  Table 2 contains estimates of the average effect of a regressor on expected export volume, conditional on non-zero trade, over the distribution of regressors and fixed effects for trading country pairs. That is, we calculate (for n 1 = i j =i 1{Y ij > 0}) for the continuous regressor log distance and for binary regressors, where X ij sets the binary regressor of interest to one for all (i, j) and X (0) ij sets it to zero (leaving other regressors unchanged). Again, the jackknife bias correction has an important impact on two of the effects; for example, the effect of a trade agreement on expected export volumes decreases by about a quarter (more than a one standard error change in magnitude). Note that, as is the case here, a small bias in the coefficient on some variable does not necessarily imply low bias in the corresponding marginal effect.

Testing for strategic interactions
To demonstrate the use of the jackknife for specification testing, we implement the test in (15) in a binary model for the probability of country i exporting to country j. We test for two types of strategic interaction: reciprocity, using S ij = Y ji ; and transitivity, using In both cases we construct the statistics by estimating p ij using a probit model and jackknife the statistic. Standard errors are computed using the expressions in Theorem 2. Table 3 presents the values of the statistics as well as the standardized values For the reciprocity statistic, the jackknife bias correction appears to have little effect. The statistic rejects the null of no reciprocity strongly suggesting that the existence of an export relationship from i to j increases the likelihood of i also importing from j. This is perhaps not surprising. It is important to note that the model considered in this paper allows for reciprocity so that this conclusion has no impact on the validity of model estimates. The presence of reciprocity does suggest that standard errors should be clustered at the dyad level to account for correlation between the outcomes Y ij and Y ji . In contrast, the jackknife bias correction appears to have an important effect on the transitivity statistic. The uncorrected statistic leads to a rejection of the null, and the conclusion that indirect paths of trade (exports paths from i to j through a third-party country) are associated with a lower probability of a direct export relationship than is expected given the model. However, the jackknifed statistic is close to zero so that we do not reject the null hypothesis that trade decisions are bilateral in nature.

Comparison with conditional logit estimator
As a comparison with an existing approach to the incidental parameters problem in networks, we consider estimating a logit model for the probability that country i exports to country j. Under the logit specification, it is possible to estimate the common parameters in the model by first removing the fixed effects. This approach has been suggested by Graham (2017) for an undirected network, and Jochmans (2018) for a directed network. The conditional logit estimator works by forming difference-in-differences style contrasts among sets of four nodes (tetrads). Let Given the logistic specification P (Y ij = 1|X, β, α, γ) = exp(β X ij + α i + γ j ) 1 + exp(β X ij + α i + γ j ) and conditional independence of outcomes across dyads, we may write P (z ij,kl = 1|z ij,kl ∈ {−1, 1}, X, β, α, γ) = exp(β r ij ) 1 + exp(β r ij ) That is, conditional on the event z ij,kl ∈ {−1, 1}, the outcome z ij,kl follows a logit model without any fixed effect parameters. This allows us to estimate the common parameter β by a standard logit regression of z ij,kl on r ij,kl in the subsample of z ij,kl ∈ {−1, 1}. Estimates from this model, as well as the jackknife bias-corrected estimates, are shown in Table 4.
Interestingly, although the jackknife bias correction suggests little bias in the original logit parameter estimates, the conditional logit estimates are significantly different. For example, while the coefficient on distance is similar across the MLE and jackknife estimates (-1.34 and -1.31 respectively), the conditional logit estimate differs by more than three standard errors.
There may be a number of reasons for such a discrepancy between the estimates. One possibility is that the model is misspecified in the sense that the correct link function is not the logistic function. In this case, maximum likelihood estimates a set of pseudo parameters that represent the parameters that minimize the Kullback-Leibler distance between the logit model and the true model (White, 1982). Since the jackknife theory does not rely on any information equalities, the jackknife bias correction is asymptotically unbiased for these pseudo parameters. In contrast, the conditional logit estimator estimates a different set of pseudo parameters, those that minimize the KL distance conditional on z ijkl ∈ {−1, 1}, for a logit regression of z ijkl on r ijkl . To the best of our knowledge, there is no reason to suspect that these two sets of pseudo parameters would coincide.
An alternative explanation is that the assumption of independence between dyads is violated in the data. In this case, the likelihood function for the conditional logit will be incorrect, since the identity (22) will no longer hold. In this setting the bias correction given by the jackknife estimator is also likely to be incorrect since it will not account for bias terms generated by the dependence across dyads.

Simulations
Here I demonstrate the effectiveness of the jackknife in simulations. I repeat the simulation design of Dzemski (2019), which has also been used in a number of other network papers. The binary outcome Y ij is determined by where θ = 1 and ε ij ∼ N (0, 1). Individual i is characterized by the binary scalar X i = 1 − 2 × 1{i is odd}, and the homophily variable is given by X ij = X i X j , i.e. it is one for pairs with the same sign and minus one for pairs with opposing signs. The fixed effects are given by which is a sequence from C l N to C u N . The value of (C l N , C u N ) is intended to control the sparsity of the network, and we consider four choices, shown in Table 5.
In the balanced setting ('bal'), fixed effects range between ± log log N , generating a dense network in which around half of all links are formed. Subsequent settings feature increasingly sparse networks in which some nodes remain well connected, while others make few links. In the sparsest setting ('log') only around 3 per cent of all links are formed. I simulate the model 1000 times and compute the MLE, analytical bias-corrected estimate, and both the standard and weighted jackknife bias-corrected estimates for each fixed effect distribution. Table 6 presents the bias, standard deviation, and rejection rates for a 5 per cent test of the null hypothesis β 0 = 1 for each estimator. As expected, the MLE is biased, with the size of the bias increasing in the sparsity of the network. In each case, the bias is around one standard deviation in size, resulting in substantial over-rejection. In contrast, the jackknife estimator is approximately unbiased in the first two designs. In the sparsest design the jackknife estimator does not appear to perform well, and shows similar bias to the MLE. The weighted jackknife performs well, even in the sparse designs; it is unbiased in the first    three designs, and removes the majority of the bias even in the sparsest design. In each case, the jackknife estimators have smaller standard deviation than the MLE (with the exception of the standard jackknife in the sparsest design). Rejection rates are at or below the nominal level in three of the four settings, although in the sparsest setting all estimators over reject.

Weighted jackknife
To help explain the properties of the weighted jackknife, we investigate how the weight given to each leave-out estimate β (k) correlates with the contribution of that estimate to the total error of the jackknife estimate. Define as the contribution to the error β J − β 0 from a single leave-out estimate β (k) . The total error of the standard jackknife is simply the average of these errors β J − β 0 = 1 N −1 k e (k) , while the weighted jackknife gives differing weights to different leave-out estimates. Figure 3 plots

versus estimates ( β (k) ) in leave-out samples
The y-axis is which measures the absolute error of the jackknife estimator using the single leave-out estimate. The x-axis is the weight given to that leave-out estimate in the weighted jackknife.
the absolute value of these errors against their weights in the weighted jackknife across the four simulation designs.
For the balanced design (upper left panel), which is the densest network, the weights are concentrated around 1 N −1 ≈ 0.02 with very little variation. In this case, the weighted jackknife is almost identical to the standard jackknife, as is clear in Tables 6 and 7. As the level of sparsity in the network increases, so does the dispersion in weights and estimates across the leave-out samples. In the sparsest design, C l N = − log N (lower right panel), the weights vary considerably, with many close to zero. Note that the sparser designs also exhibit large outliers, with extremely large errors, but that these outliers typically receive weights close to zero. It is this feature that appears to drive the success of the weighted jackknife over the standard jackknife in the sparser designs (C l N = −{log 1/2 N, log N } ) in Tables 6 and 7.
Simulations for leave-l-out style jackknife Table 8 reports the results of simulations of the leave-l-out style jackknife in a network of N = 101 nodes. The model is the same as that in the previous section, with fixed effect distributions given in Table 5. Results are shown for l = {5, 10, 20, 50}, which corresponds to 20, 10, 5, and 2 leave-out estimates in each case. The leave-l-out jackknife performs well,  even for l as large as 50 in the dense network scenarios. In the sparsest scenario, the jackknife performs reasonably well for l as large as 20, but appears to break down for l = 50. For large l = 50, half of all observations in the network are dropped for each leave-out estimation, which appears to create some issues when the network is sparse.

Comparison with other jackknife bias corrections
Finally, we compare the jackknife suggested in this paper to previous suggestions. Specifically, Cruz-Gonzalez et al. (2017) suggest two jackknife bias corrections for network models of the type considered in this paper. 9 The first bias correction is based on a split-sample approach. Divide the agents into two halves, A 1 and A 2 . Define β α,γ/2 = 1 2 ( β α,γ∈A 1 + β α,γ∈A 2 ), where β α,γ∈A 1 is the estimator that uses only observations in which the receiver is in the first set of agents A 1 , and β α,γ∈A 2 uses only observations in which the receiver is in A 2 . Similarly, define β α/2,γ = 1 2 ( β α∈A 1 ,γ + β α∈A 2 ,γ ) as the average of the two estimators that split the sample based on sending agents. A split-sample jackknife is given by The second bias correction is based on dropping all observations associated with a particular 9 Cruz-Gonzalez et al. (2017) suggest these jackknife corrections for the network model although do not prove their validity. Establishing the validity of the jackknife corrections requires the higher-order asymptotic expansions that are derived in this paper.
Both jackknife corrections differ from β J in (8), only β J preserves the distribution of fixed effects in each leave-out estimate. This appears to be an important property for the jackknife to perform well in settings with substantial unobserved heterogeneity. Table 9 reports results from simulations of the same model discussed above for the standard and weighted jackknife estimators as well as the split-sample and double corrections suggested by Cruz-Gonzalez et al. (2017). As is clear from the results, although the split-sample and double corrections appear to work well in the densest network DGP, they do not perform well in settings with more heterogeneity and sparser networks. In particular, the split-sample correction has much larger variance, and removes less bias than the leave-one-out style corrections (Hughes and Hahn (2020) derive higher-order bias and variance expressions that explain this phenomenon in the panel setting).

Conclusion
This paper presents a new method for bias correcting nonlinear dyadic network models with fixed effects. We provide a novel formulation of the jackknife method that applies to networks with both sender and receiver fixed effects. The jackknife method provides an 'off-the-shelf' procedure for bias correction that is easy to apply, and applicable to a wide set of models.
It allows for discrete multivalued and continuous outcome variables, and is able to obtain estimates of average effects and counterfactual outcomes.
In addition, we show how the jackknife can be used to bias correct averages of functions that depend on multiple observations, including dyads, triads, and tetrads in the network. These averages can be used to produce a wide array of test statistics for the presence of strategic interactions in the network, such as reciprocity or transitivity. In simulations, we show that the jackknife performs well, even in relatively low density networks, and outperforms previous suggestions for jackknife procedures.
There are a number of interesting areas in which the work in this paper might be usefully extended. The jackknife procedure developed in this paper applies to data on a single observation of a network. It would be interesting to extend these results to networks observed over multiple time periods, perhaps with the addition of time fixed effects. It is expected that a similar jackknife procedure could also be used in this setting, with appropriate splitting across the time dimension to account for dynamics in the network. The jackknife procedures proposed in this paper may also be useful in the interactive fixed effect model of Chen et al. (2021), and establishing their validity in this setting would also be useful.

A Notation and norms
The notation in the appendices follows that in the main paper. That is, denote partial derivatives of the objective function using subscripts, so that ∂ β L(β, φ) denotes ∂L(β, φ)/∂β and so on, where φ = (α , γ ). When functions are evaluated at β 0 , φ 0 the dependence on these arguments is dropped. We also use ij as shorthand for the function (Z ij , ·), with be the first derivatives of the objective function. We also write for the negative of the Hessian with respect to the fixed effects, a 2N × 2N matrix.
We follow Fernández-Val and Weidner (2016) (FVW) in using the Euclidean norm · for dim β vectors, and the norm induced by the Euclidean norm for matrices and tensors, i.e.
In the proofs we sometimes take β to be a scalar to simplify notation, although the results apply to any vector of fixed size. Since the number of fixed effect parameters in the model grows with N , the choice of norm for dim φ vectors and matrices is important. Following FVW, we choose the q -norm for dim φ vectors and the corresponding induced norms for matrices and tensors See FVW for more details on these norms.
We define the sets B(r, B Asymptotic expansions B.1 Verifying Assumption B.1 in Fernández-Val and Weidner (2016) The results in the paper are based on asymptotic expansions of the objective function. Fernández-Val and Weidner (2016) (FVW) derive expansions for a general class of Mestimators with multiple incidental parameters, which includes the model studied here. In order to determine the properties of the jackknife estimator, I extend the expansions in FVW to higher-order, which requires additional conditions on the number of moments and derivatives of the objective function that exist. Derivations of higher-order terms and their bounds are quite long and largely similar to the derivations in FVW, and so are provided in the Supplementary Appendix. This appendix contains analysis based on the first-order expansions and focuses on results related to the jackknife, which are of most interest. To do so, I begin by verifying the conditions in Assumption B.1. of FVW, which will allow us to use some important results in that paper, including consistency of the common parameter and vector of fixed effects.
As is the paper, we useĒ to denote expectation conditional on exogenous covariates and fixed effects.
For part (iv), by (23), the condition holds by Assumption 1 (v). Condition (v) follows similarly to the proofs in FVW for the panel case. For example, we have by Assumption 1 and similarly for ∂ β k γ L q , and hence ∂ β k φ L q . For part (vi) the proofs in the panel case again carry over, for example see the proof in Lemma A.2 of Dzemski (2019). Inspection of that proof shows that we can get the bound H = O p (N 2 − 1 2 ), for > 1 32 , which will be useful in the higher-order approximations.
Given Assumptions B.1 (i), (ii), (iv), (v), and (vi) we can apply Theorem B.3 in FVW to establish consistency of the estimates for common parameters and fixed effects, and to establish the bounds

B.2 Asymptotic expansions
The next lemma gives the asymptotic expansion for the estimated fixed effects and common parameters. It is in part a restatement of Theorem B.1 in Fernández-Val and Weidner (2016), but the remainder terms are split into two parts. The expressions for the remainder terms are given in the Supplementary Appendix, and are used to derive the properties of the jackknife bias correction. The proof of Lemma 1 is provided in the Supplementary Appendix. Lemma 1. Let Assumption 1 hold. Then,

C Jackknife results for β
Here we allow for a more general construction of the leave-out sets, but impose two important conditions.
Condition 1. Let I k for k = 1, . . . N k be a partition of the observations in a network of size N . Define 1 k ij = 1{(i, j) ∈ I k } as an indicator that the observation (i, j) is not included in the k-th set. We impose the following conditions on the sets: Condition 1 imposes two important constrains on the sets I k : (i) that they are mutually exclusive, such that every edge (i, j) appears in exactly one of the sets, and (ii) that each set contains exactly one observation related to each of the N sender fixed effects α and each of the N receiver fixed effects γ. The first condition ensures that all observations are used equally in the jackknife, and is important for showing that the asymptotic variance of the estimator is not affected by the jackknife. The second condition ensures that each fixed effect parameter is affected equally in the leave-out sets, and thatH (k) is well-defined and positive definite. We assume in the proofs below that N k = N − 1 and that we are using the leave-one-out jackknife. All of the results still hold for the leave-l-out style jackknife for fixed l.
Some additional notation will also be useful for studying the jackknife estimates. Let A be some statistic and A (k) the same statistic in the k-th leave-out sample. We define the jackknife operator J as Additionally, we define a set of indicators that count the number of unique leave-out sets I k that a group of edges (i 1 , j 1 ), . . . , (i t , j t ) are contained in. Let Assume that the researcher randomizes the ordering of the N node labels by choosing, with equal probability a labelling from the N ! possible orderings. Since the ordering of nodes is random, the indicator 1 k ij is a random variable, with mean equal to the probability that (i, j) ∈ I k across each of the possible orderings. There are N × (N − 2)! ways to order the N nodes keeping (i, j) ∈ I k for any fixed k, so that E[1 k ij ] = N −2 N −1 . The following lemma states that the expectation, over the randomness induced by the ordering of nodes, of sums in the leave-out sets is equal to that in the full sample.
Lemma 2. Let 1 k ij satisfy Condition 1, and define the sums over random variable A ij in the full-sample and k-th leave-out sample Proof. The proof follows simply from the fact that E[1 k ij ] = N −2 N −1 , where the expectation is taken over the randomness induced by the random ordering of nodes. Note that the ordering of nodes is independent of any other randomness in the sample, so that Then we haveĒ as required.

C.1 Jackknifing first-order expansion terms
We next use this result to show the impact of the jackknife operator on the first-order expansion of N ( β − β). From the expansion in Lemma 1, we have that a first-order approximation is given by where The next lemmas demonstrate the effect of the jackknife on general sums of the forms in U (0) and U (1) .
Lemma 3. Let 1 k ij satisfy Condition 1. For A ij a mean-zero random variable, let Proof. First note that by Condition 1, k 1 k is = N − 2. Then, k A (k),i = s =i A is so that Lemma 4. Let 1 k ij satisfy Condition 1. For A ij a mean-zero random variable with bounded fourth moment, let and let B and B k be similarly defined vectors of sums involving mean-zero random variables where M is a non-random matrix that has O p (1) elements on its diagonal and O p (N −1 ) off-diagonal terms. Then we have: Proof. The most common choice of M will beH −1 , which satisfies the conditions for M by Assumption 1 and the approximation property in (23). We show the proof usingH −1 , but note that it holds for any M satisfying the conditions stated above. We have Recall that I 1 (ij)(st) is one whenever (i, j) and (s, t) are contained in the same I k , and so Then, we have Similar computations for the other three elements gives We can decompose J 0 as Then we havē Where j =i ∂ π ij . Note that the jackknife operator is linear so that, since (∂ βφ L )H −1 are fixed across leave-out samples (by Lemma 2), Now, for U (1) we can appeal to Lemma 4. For the first term, (∂ βφ L )H −1 S we set A ij = ∂ βπ˜ ij and B ij = ∂ π ij . For the second term, (∂ βφ L )H −1HH−1 S we note that and similarly for elements of (∂ βφ L )H −1H ·,γ . So we can let and B ij = ∂ π ij in Lemma 4. Note that we have A ij and B ij mean zero and the moment condition also holds by assumption. For the final term, we begin by demonstrating that bothH −1 g (∂ βφ φgL )H −1 andH −1 g (∂ φφ φgL ) H −1 (∂ βφL ) gH −1 are matrices that satisfy the requirements for M in Lemma 4. The proof for the second term is shown, with the result for the first term following nearly identically. Firstly, Taking the first of these terms, the (i, j) element is given by Identical results apply to the other elements in g (∂ φφ φgL ) H −1 (∂ βφL ) g and hence we can conclude that the matrix has O p (1) diagonal elements and O p (N −1 ) off-diagonal elements. It then follows that the same is true ofH −1 g (∂ φφ φgL ) H −1 (∂ βφL ) gH −1 . Then, we can apply Lemma 4 with A = B = S to give the result.
Lemma 6. Let Assumption 1 hold, and let β J be the either the jackknife, leave-l-out jackknife, or weighted jackknife estimator. Then, a first-order approximation to the estimator is given byW Proof. Recall from Lemma 1 that SinceW N is fixed across leave-out samples (Lemma 2), we can focus on the jackknife operator applied to the RHS. By Lemma 5 we have that J[ since each remainder term in the above is o p (N −1 ).

C.2 Approximation of W N
The next two results show that the sample version of the Hessian for the common parameters β is consistent, and that it is approximately the same across leave-out samples.
Lemma 7. Let Assumption 1 hold. Then, for ≥ 1 32 The first term is 1 For the remaining term, we can decompose it as For the leave-out term note that the moment bounds in Assumption 1 (iii) imply identical bounds in the leave-out samples, simply by replacing Lemma 8. Let Assumption 1 hold, then for all k Proof. We prove the first statement, since the proof of the second is identical. A first-order Taylor expansion of W N gives whereβ andφ are intermediate values between (β 0 , φ 0 ) and ( β, φ). From Assumption 1 and the bounds in Assumption B.1 of FVW, we have that W N is differentiable with derivatives that are O p (1), since (1) follows simply by expansion of the matrix product and the properties of H −1 in (23). Similarly, Then, since β − β 0 → 0 and φ − φ ∞ → 0 by (24), we get the result.

C.3 Proof of Theorem 1
From Lemma 6 we have that The result then follows from a standard CLT argument, noting that D β ij + D β ji are independent over (i, j). For the weighted jackknife, we use the fact that Lemmas 7 and 8, and the triangle inequality give W J −W N = o p (1) and hence and hence the weighted jackknife is equal to the standard jackknife to first-order.
To show consistency of the plug-in estimator Ω N , we note that by Assumption 1 (iii), D β ij has a first-order Taylor approximation that is a continuously differentiable function of the parameters. Then by the continuous mapping theorem and the consistency results β−β 0 → D Jackknife results for average effects We begin by stating a first-order asymptotic expansion for the average effect estimator that will be used in the proof of Theorem 2. The proof of this result is provided in the Supplementary Appendix. Lemma 9. Let Assumptions 1 and 2 hold. Then In order to establish an equivalent asymptotic expansion for the leave-out estimators ∆ (k) we need to determine the value of expectations in the leave-out samples. The next lemma does this, and is analogous to Lemma 2, which states the same result for averages over single observations.
Lemma 10. Let 1 k ij satisfy Condition 1, and define 1 k λ = (i,j)∈λ 1 k ij for λ a set of r observations (i, j). Then, for sums Proof. We will prove that E[1 k λ ] = N −r−1 N −1 + O(N −2 ) from which the statement in the lemma follows sinceĒ Let I λ be equal to the number of leave-out sets I k spanned by the observations in λ. Since any observation is equally likely to appear in any leave-out set, We begin by arguing that P (I λ = r) = 1 − O(N −1 ), that is, the sets λ do not contain multiple observations from within the same leave-out set with probability approaching one. It is sufficient to consider sets λ that contain p = 2r unique agents, since this case represents the most likely scenario for λ containing observations in the same leave-out set as within a leave-out set no two senders can be the same, and no two receivers may be the same. There are N ! (N −2r)! possible choices for the ordered set of agents in λ. Then, there are at most N (N − 1)(N − 2)(N − 4) · · · (N − 2r + 2)(N − 3r + 2) orderings of agents in λ which span r different leave-out sets (N (N − 1) choices for the first observation (i, j), then (N − 2)(N − 4) choices for the second observation (s, t) since (s, t) cannot belong in the same leave-out set as (i, j), and so on). This gives which is the product of 2r ratios each of which is equal to 1 − O(N −1 ), which implies P (I r λ = 1) = 1 − O(N −1 ).
Next, note that whenever I λ = s, we have 1 k λ = 1 only if I k is not one of the s leave-out sets spanned by λ. This happens with probability N −2 and so Lemma 11. Let λ be a set of r observations (i, j) involving p unique agents, and Λ N be the collection of all such λ formed by permuting the agents in λ. Then, under Assumption 2, Proof. Let λ α denote the set of p α sender agents in the observations within λ, and λ γ the set of p γ receiving agents. There are |Λ N | = N ! (N −p)! ways of selecting the p agents in λ. Among these permutations, agent i is a sender p α (N −1)! (N −p)! times, while node j is the receiver p γ (N −1)! (N −p)! times. Using this, the first derivatives of∆ N with respect to the fixed effects are where the O p (N −1 ) statements come from the fact that p α N . An identical result applies to the diagonal elements of ∂ φφ ∆ N , i.e. ∂ α i α i∆ N = O p (N −1 ) and ∂ γ j γ j∆ N = O p (N −1 ) since they are sums over the same sets of λ. Also, if the presence of i as a sender agent implies that i is also a receiver in λ (e.g. the cyclic triangle {(i, j), (j, k), (k, i)}) then it will be the case that ∂ α i γ i∆ N = O p (N −1 ) also (if this is not true it will be lower order).
Next, consider the off-diagonal components of ∂ φφ ∆ N . If p α = 1 then ∂ α i α j∆ N = 0, otherwise, there are pα 2 (N −2)! (N −p)! permutations that contain both i and j as senders. Similarly, for p γ ≥ 2, there are pγ 2 (N −2)! (N −p)! permutations that contain both i and j as receivers. Finally, there are at most p α p γ (N −2)! (N −p)! permutations in which i is a sender and j a receiver (this is an upper bound since with i in a particular sender position, not all receiver positions may be valid for j). This, along with Assumption 2, gives the results which demonstrates the lemma.
Lemma 12. Let 1 k ij satisfy Condition 1 and let 1 k λ = (i,j)∈λ 1 k ij . Let A ij be a mean-zero random variable with bounded fourth moment, and define and let B and B k be defined as for mean zero B λ with bounded fourth moment. Assume that A ij is independent of A st for (i, j) ∈ {(s, t), (t, s)}, and independent of B λ whenever λ does not contain either (i, j) or (j, i). Define the jackknifed term where M is a non-random matrix that has O p (1) elements on its diagonal and O p (N −1 ) off-diagonal terms. Then we have: Proof. The most common choice of M will beH −1 , which satisfies the conditions for M by Assumption 1 and (23). We show the proof usingH −1 , but note that it holds for any M satisfying the conditions stated above. We have Let Λ is = {λ : (i, s) ∈ λ} be the set of λ containing observation (i, s). The full sample and leave-out versions of the first term are Let I λ be equal to the number of leave-out sets I k spanned by the r observations in λ so that k 1 k λ = N − I λ − 1. As shown in the proof of Lemma 10, |{λ : I λ < r}|/|Λ N | → 0, that is, the fraction of sets λ that contain two or more observations in the same I k is a vanishingly small. Using this, we have

similar computations for the other three elements gives
Next, recall that Γ ijst = O p (1) whenever i = j or s = t, and is O p (N −1 ) otherwise. Note thatĒ[J 0,3 ] = 0 since λ does not contain (i, s) or (s, i). Taking expectations we get , and so on applying the results on the size of sets Λ ij , I λ < r, and Γ ijst .
Then, An almost identical analysis applies to J 0,1 and J 0,2 , giving the result J 0 = o p (1).
The next lemma states the first-order approximation for the jackknife bias-corrected average effect estimator.
Lemma 13. Let Assumptions 1 and 2 hold and let ∆ J be the jackknife bias-corrected estimator in (12). Then, Proof. We can write an expansion for the leave-out estimate where R (k),∆ = o p (1) is the version of R ∆ in the leave-out sample andR (k),∆ = o p (N −1 ) is the leave-out version ofR ∆ combined with the error from replacing terms like ∂ β∆(k) with ∂ β∆N (i.e. applying the result in Lemma 10). Using the expansion for the leave-out estimate, we can apply the jackknife operator to each line above.
For the first term, we apply the result in Lemma 5 to give Similarly, Lemma 3 implies that jackknifing the second term gives N (∂ φ ∆ N )H −1 S.
For the third term, we note that by Lemma 11 we can apply Lemma 4 with M =H −1 , A = S, and B = N (∂ φ ∆ N )H −1H , and apply Lemma 12 with M =H −1 , A = S, and either B = N (∂ φ ∆ N ).
For the fourth term, we first show that satisfies the conditions for M in Lemma 4, from which we will be able to conclude that the jackknifed term will be o p (1). The above expression is non-random (conditional on exogenous regressors and fixed effects) and so we must demonstrate that it is a 2N × 2N matrix with O p (1) diagonal elements, and O p (N −1 ) off-diagonal elements. Note that if two 2N × 2N matrices both have O p (1) diagonal elements and O p (N −1 ) off-diagonal elements, then their product also shares this property. Since this is true ofH −1 (see Lemma D.1 in Fernández-Val and Weidner (2016)), it remains to demonstrate this fact for the terms N (∂ φφ ∆ N ) and N g (∂ φφ φgL )[(∂ φ ∆ N )H −1 ] g .
By Lemma 11 we have that N (∂ φφ ∆ N ) has O p (1) diagonal elements and O p (N −1 ) offdiagonal, with the possible exception of the elements ∂ α i γ i∆ N . However, this still implies that NH −1 (∂ φφ ∆ N )H −1 satisfies the condition. For N g (∂ φφ φgL )[(∂ φ ∆ N )H −1 ] g , diagonal elements are given by (for i ≤ N ) which is O p (1) since by Lemma 11 and Lemma D.1 in Fernández-Val and Weidner (2016), and similarly for i > N . Off-diagonal components can be shown similarly, e.g. for i < N and j > N we have which is O p (N −1 ). Finally, in the Supplementary Appendix (S.4) it is shown that J[R ∆ ] = o p (1), and byR ∆,(k) = o p (N −1 ) for each k (and in the full sample) we have that (N − 2)R ∆,(k) = o p (1) and (N − 1)R ∆ = o p (1) so that the jackknifed version of this term is also o p (1).

D.1 Proof of Theorem 2
We can decompose ∆ J −∆ N into From Lemma 13 we have Some tedious matrix algebra shows that this expression is equivalent to where θ = (β, φ ), andh ij = −N (∂ θ∆N )(∂ θθ L ) −1 ∂ θ ij . Next, lets ij = 1 which is asymptotically normal by a standard CLT. We now determine the variance of this term.

V ar
To compute this, first note that s t =sĒ [h ijhst ] =Ē[h ij (h ij +h ji )]. Also, where s ij = (N −p)! (N −2)! λ∈(Λ ij ∪Λ ji ) (m λ −m λ ). Let D(λ) be the set of dyads formed from the observations in λ, i.e. if (i, j) ∈ λ then (i, j), (j, i) ∈ D(λ). Assume that λ and λ both contain the observations (i, j),(i, k). There are O(N p−3 ) sets λ containing the corresponding dyads, so that there are O(N 2p−3 ) such λ, λ pairs (O(N 3 ) choices of (i, j, k) and O(N p−3 ) choices for each of λ and λ ). Similarly, there are O(N 2p−4 ) λ, λ pairs that share two dyads made up for four agents, (i, j),(k, l). for h ij =h ij +h ji . The asymptotic variance of N ( ∆ J −∆ N ) is given by the limit of this expression. Assumptions 1 (iii) and 2 (ii) guarantee that both h ij and s ij have first-order approximations that are continuously differentiable in the parameters, so that the continuous mapping theorem and the consistency results β − β 0 → 0 and φ − φ ∞ → 0 in (24), imply consistency of the plug-in estimator for V ∆ (see for example Lemma S.1 in FVW).

D.2 Proof of Theorem 3
We begin with a U-statistic representation of∆ N , which will allow us to apply standard asymptotic results on U-statistics. We have defined m to be a function of the sets λ, which depend on an ordered set of p agents. For example the transitive triangle λ = {(i, j), (i, k), (k, j)} depends on the agents {i, j, k} in a non-symmetric manner. We first rewrite∆ N to be a sum over functions that are symmetric in agents. Denote the set of agents in λ by N (λ), and let η = {i 1 , . . . , i p } be some set of p agents. Then we may define τ η = {λ : N (λ) = η} as the collection of all λ that contain the same set of agents. Note that |τ η | = p!. We have, form =m − E[m] where u τ = 1 p! λ∈τm λ . The variable u τ is symmetric function of {β, X i , α i , γ i } for p agents i. For example, there are 3! = 6 possible transitive triangles using agents {i, j, k} so that u is the average of the function m evaluated at these 6 different triangles. Assuming that the {X i , α i , γ i } are i.i.d. over agents,∆ N − δ is a U-statistic of order p and we apply standard theory on such statistics to compute its asymptotic distribution. As in Theorem 12.3 in van der Vaart (1998) we have where, for τ and τ sharing exactly one agent in common, The variance of t i is given by ways to choose the remaining p − 1 agents in τ so that τ and τ share only agent i in common. The first term is therefore O(1). Now assume τ and τ share two agents in common (one of which is i). There are N − 1 choices for the second common agent, N −2 p−2 ways to choose τ containing i and the second common agent, and N −p p−2 ways to choose the remaining agents in τ . This implies that the sum for τ and τ with two agents in common is O (N −1 ). Similarly, the sums for three agents in common are O(N −2 ) and so on.
This implies that the variance of t i converges to ζ 1 . We can alternatively express this An estimator of p 2 ζ 2 1 is therefore given by and m λ is a plug-in estimator form λ . Assumption 2 and consistency of parameters φ − φ 0 ∞ → 0, β − β 0 → 0 ensures consistency of m λ and hence consistency of V δ .