Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
120,430 characters · 20 sections · 102 citation commands
Estimation of Peer Effects in Endogenous Social Networks: Control Function Approach
The ways in which interconnected individuals influence each other are usually referred to as peer effects. One of the first to formally model peer effects is Manski1993a. He proposes the linear-in-means model, in which an individual's action depends on the average action of other individuals and possibly also on their average characteristics. Manski1993a assumes that all individuals within a given group are connected. Later literature allows for more complex patterns of connections, in which an individual might be directly influenced by a subset of the group. Examples are Bramoulle2009, Lee2010, Lee2007 among others. Models of peer effects have been applied in various areas, such as education, health and development. Examples of applications are found in recent review papers such as Blume2010, Manski2000, Epple2011, brock2001interactions and Graham2011.
Many models considered in earlier literature assume that connections between individuals are independent of unobserved individual characteristics that influence outcomes. However, assuming exogeneity of the network or peer group is restrictive in many applications. For example, consider the following widely studied empirical application of peer effects: peer influence on scholarly achievement. The assumption that friendships are exogenous in the outcome equation for scholarly achievement means that there are no unobserved variables that influence both friendship formation and individual grades. However, even if a study controls for observable individual characteristics such as gender, age, race and parents' education, it is likely to omit factors that influence both students' choice of friends and their GPA; for example parental expectations, psychological disorders, or non-reported substance use. For more examples of endogenous peer groups see brock2001interactions, Weinberg2007, Shalizi2012 and hsieh2016social, among others.
In this paper we propose a method for estimating a linear-in-means model of peer effects, where the peer group is defined by a network that is endogenous in the outcome equation. Our model allows for correlation between the unobserved individual heterogeneity that impacts network formation and the unobserved characteristics of the outcome. For this, we use a dyadic network formation model that allows the unobserved individual attributes of two different agents to influence link formation, and in which links are pairwise independent conditional on the observed and unobserved individual attributes. The network formation we consider in the paper is dense and nonparametric.
The main contributions of the paper are methodological. First, given the endogenous peer group formation, we show that we can identify the peer effects by controlling the unobserved individual heterogeneity of the network formation equation. Second, we propose an empirically tractable implementation of the control function, whose functional form is not parametrically specified. For this, we propose two approaches, one based on an estimator of the unobserved individual heterogeneity and the other one based on the average node degrees of the network.\footnote{We acknowledge that this approach is developed based on an idea provided by one of the referees. We thank the referee.} Our estimation method is semiparametric because we do not restrict the functional form of the control function. Finally, we derive the limiting distributions of the estimators within a large single network. The main challenge of the asymptotics is handling the strong dependence of observables caused by the dense network. Other peer effects papers that have considered endogenously formed peer groups and have controlled the endogeneity via various control functions include GoldsmithP2013, hsieh2016social, Qu2015, Arduini2015a and Auerbach2016. We provide more detail on these papers in Section (ref).
The remainder of the paper is organized as follows. In Section (ref) we present a high level description of our approach and provide intuition as to its empirical applications. In Section (ref) we formally present our model. In Section (ref) we show {\color{black} how to identify peer effects using control functions.} Estimation is discussed in Section (ref), and in Section (ref) we discuss the limiting distribution of the estimator and propose standard errors. In Section (ref) we present results of Monte Carlo simulations. There we compare the finite sample performance of our two semiparametric estimators against an estimator that assumes unobserved characteristics enter in a linear way, as well as an instrumental variables (IV) estimator that does not control for network endogeneity. {\color{black} We investigate both high degree and low degree networks.} Section (ref) concludes.
A word on notation: in what follows we denote scalars by lowercase letters, vectors by lowercase bold letters, and matrices by uppercase bold letters.
In this section we introduce a simple model in order to illustrate the main points of our approach. A more general model and detailed discussion of the model will follow later.
A simple peer effect model for the purpose of illustration of the main idea is
where $x_i$ is a measure of observable characteristics of individual $i$ and $d_{ij}$ is an indicator of individual $i$'s peer, so $d_{ij}=1$ if $i$ and $j$ are directly linked and $0$ otherwise. In ((ref)), the regressor of interest is the average of the characteristics of those individuals who are linked with $i$, $\frac{\sum_{j \neq i}d_{ij}x_j}{\sum_{j\neq i}d_{ij}}$. For simplicity, we assume that $x_i$ is exogenous with respect to all the unobserved components of the model; this will be relaxed later.
For the link formation, we consider the following dyadic network formation model,
where $a_i$ and $a_j$ are unobserved individual specific characteristics, $u_{ij}$ is a link specific component, and $g(\cdot,\cdot)$ is some function. It should be noted that this model of network formation does not allow for network effects in link formation, as a link between $i$ and $j$ only depends on the characteristics of $i$ and $j$.
The unobserved individual characteristic $a_i$ can be interpreted as social capital that increases the likelihood of forming a link. Depending on the context this could be factors like trustworthiness, socioeconomic status, or outspokenness.
For example, fafchamps2011 measure the risk sharing links between households in Tanzania and they construct links between households based on the question whom individuals could “personally rely on for help.” fafchamps2007risk examine the formation of risk-sharing networks using data from the rural Philippines. Banerjee2013 examine how participation in micro-finance diffuses through a social network which they measure using lending and trust. In these settings, we can think of $a_i$ as a measure of individual trustworthiness and integrity in financial matters. Ductor2011 analyze whether knowledge of a researcher's co-authorship network is helpful in predicting his or her productivity. In this setting $a_i$ can be interpreted as some unobserved productivity trait that induces the researcher to have more coauthors, and also to be more productive at writing papers.
The key feature of the peer effect model ((ref)) and ((ref)) is that individual $i$'s unobserved characteristic $a_i$, which impacts link formation, is correlated with $v_i$, $i$'s unobserved characteristic that affects the outcome $y_i$. For example, $a_i$ could be an unobserved component that affects a researcher's publication rate $y_i$, and also his or her co-authorship relationships, $d_{ij}$. Alternatively, we can think of a situation where there are two types of agents: popular and unpopular. The popular agents are more likely to be friends with other agents, and popular agents have better outcomes even in the absence of a peer effect. Then the peer formation $d_{ij}$ becomes correlated with the unobserved component $v_i$ of the outcome, and, as a consequence, the regressor of the peer effect, $\frac{\sum_{j \neq i}d_{ij}x_j}{\sum_{j\neq i}d_{ij}}$, becomes endogenous.
In this paper we use a control function method to handle the endogenous peer group problem. Let $\mathbf{D}_N$ be the $N \times N$ adjacency matrix that describes the network links $d_{ij}$. Suppose that the unobserved characteristics $(a_i,v_i)$ and $u_{ij}$ are randomly drawn over $i$ and $(i,j)$, respectively. Also assume that $u_{ij}$ is independent of $(a_i,v_i)$. Then, for any $i \neq j$, the link $d_{ij} = \mathbb{I}(g(a_i,a_j)\geq u_{ij})$ and $v_i$ are dependent only through $a_i$. Therefore, controlling for $a_i$, the network $\mathbf{D}_N$ and $v_i$ become mean independent, that is,
Suppose that we observe $a_i$. Consider the outcome equation which controls for $a_i$ nonparametrically, \[ y_i = \beta^0 \left( \frac{\sum_{j \neq i}d_{ij}x_j}{\sum_{j\neq i}d_{ij}} \right) + h(a_i) + \varepsilon_i, \] where $\varepsilon_i := v_i - h(a_i)$. Once we control the endogeneity of the network with $a_i$, then the regressor of the peer effect becomes exogenous, and we can estimate the peer effect coefficient $\beta^0$ using the conventional partially linear regression estimation method (e.g. Robinson1988).
However, in most empirical applications, $a_i$ is not observed. Then the question becomes how to implement the control function. In this paper, as the main methodological contribution, we propose the following two procedures. Both procedures are implemented with a single snapshot of an observed network.
The use of degree as a control function requires much fewer restrictions on the specification of the network. Intuitively, the unobserved node (or individual) fixed effects $a_i$ control for heterogeneous degree sequences. Therefore, from an economic point of view, what needs to be controlled is the agent's degree, which validates the control function approach that uses $\rm{deg_i}$. This approach does not require a specification of the specific structural model explaining heterogeneous degree sequences. Consistent estimation of $a_i$ usually requires a specific functional form. For example, Graham2017 assumed an additive model and ChenFernandez-ValWeidner2018 require an interactive form. However, there is a disadvantage in the degree approach. The degree approach cannot identify the coefficient of the observed exogenous regressor if the same regressor also impacts the network formation.
In Section (ref), we generalize the simple model ((ref)) by allowing for an additional peer effect, $\frac{\sum_{j \neq i}d_{ij}y_j}{\sum_{j\neq i}d_{ij}}$, known as the endogenous peer effect, which measures the effects of the outcomes of the peer group on an individual outcome. In this case we have to deal with two kinds of endogeneity in the peer effect regressors: one from the endogenous regressors $y_{j}$ and the other one from the endogenous peers $d_{ij}$. In Section (ref), we also generalize the dyadic network formation model by introducing a dyadic component based on observed individual characteristics. We provide application examples of the general model and discuss its features there. The identification of the peer effects in the general model will be discussed in Section (ref). In Section (ref) we shows how to implement the two aforementioned estimation methods in the general framework. In the appendix we provide the regularity conditions that are required for the asymptotic results of the paper. All the technical proofs and comprehensive Monte Carlo simulation results are found in the Online Supplement material which is available in JohnssonMoon2019.
Closely related papers that adopt a control function approach include GoldsmithP2013, hsieh2016social, Qu2015, Arduini2015a and Auerbach2016. Our paper adopts a frequentist approach based on a nonparametric specification of the network formation, while GoldsmithP2013 and hsieh2016social use the Bayesian method based on a full parametric specification of the network formation and the outcome equation. Like our paper, Qu2015 assume the network (spatial weights in their model) to be endogenous through unobserved individual heterogeneity. However, our paper is different from Qu2015 in many aspects. They consider sparse network formation models while we consider a dense network. They restrict the functional form of the control function to be linear, while we impose no restriction on the functional form. The two papers propose different implementations of the control function. Also, in GoldsmithP2013, unobserved components account for homophily in link formation, whereas in our setup they mainly drive degree heterogeneity but are allowed to account for homophily as well, as in the example ((ref)).
Our paper is different from Arduini2015a regarding the main source of the endogeneity of the network and the form of the control function. Arduini2015a assume that the endogeneity of the network is allowed through dependence between the outcome equation error and the idiosyncratic network formation error, like the conventional sample selection model. This model can be interpreted as meeting opportunities being correlated with unobserved ability of the agent that affects the outcome. Arduini2015a consider control functions (both parametric and semiparametric) to deal with the selection bias problem and propose a semiparametric estimator that uses a power series to approximate selectivity bias terms. Regarding asymptotics, in both Qu2015 and Arduini2015a, the asymptotics are derived using near-epoch dependence and are based on the assumption that the number of connections does not increase at the same rate as the square of the network size.
Among the aforementioned related papers, probably the one most closely related to ours is Auerbach2016. As a result, we would like to discuss the differences between the two papers in more detail. The outcome model of Auerbach2016 is a partially linear regression model where the nonparametric component is an unknown function of the unobserved network heterogeneity,
In the simple peer effect example, the exogenous peer effect corresponds to the regressor $x_i$ above. The network formation is the same as ((ref)).
To compare the identification ideas, let's assume that $a_i \sim U[-1/2,1/2]$ and $u_{ij} \sim U[0,1]$. In this case, $d_i := (d_{i1},...,d_{in})'$ and the distribution of $d_i$ of node $i$, whose characteristic is $a_i$, is fully characterized by the link formation probability profile $g(a_i, \bullet)$.
The key condition of Auerbach2016 is that $h(a_i)$ and the the link formation distribution profile $g_i(\bullet) := g(a_i,\bullet)$ be one-to-one a.s., that is, $g(a,\bullet) \neq g(a^*,\bullet)$ a.s. if and only if $h(a) \neq h(a^*)$. Then, for any distance measure between the two profiles $g_i$ and $g_j$, $d(g_i,g_j)$, it follows that $d(g_i,g_j) = 0$ if and only if $h(a_i) = h(a_j)$.
Based on this, Auerbach2016 finds that one can control the network endogeneity by pair-wise differencing\footnote{This resembles Powell1987, Heckmanetal1998, and AbadieImbens2006.} of the observations of the two individuals, $i$ and $j$, whose network formation distributions are the same, $d(g_i,g_j) = 0$, and proposes a semiparametric estimator based on matching pairs of agents with similar columns of the squared adjacency matrix.
Notice that the identification condition of Auerbach2016 is satisfied if $g(a_i,\bullet)$ and $a_i$ have a one-to-one relation. However, our second identification is based on the condition that $a_i$ and the marginal network probability, $\int g(a_i,\tau) d \tau$, have a one-to-one relation. We admit that this condition is more restrictive than the identification condition of Auerbach2016, because our restriction is a special case of his restriction. However, as mentioned in the introduction, our identification under the stronger condition allows for the omitted variable in the peer effects equation to be nonparametrically directly estimated, which results in the peer effect estimator having the parametric convergence rate ($\sqrt{N}$). This feature is not necessarily guaranteed in the framework of Auerbach2016.\footnote{We thank one of the referees for suggesting the comparisons.}
In this section, we introduce a general linear-in-means peer effect model that extends the simple illustrative outcome model with a peer effect in ((ref)) and the simple dyadic network formation model in ((ref)).
As in Section (ref), $d_{ij}$ are the observed binary variables that measure undirected links among individuals $i\in \{1,2,\ldots,N\}$. We assume that individual outcomes are given by the linear-in-means model of peer effects
where $\mathbf{x}_{1i}$ are observed individual characteristics that affect the outcome $y_i$, $v_i$ are unobserved individual characteristics, and \[ g_{ij} = \left\{
\right. \] is the weight of the peer effects. Using the terminology of Manski1993a, $\beta_1^0$ captures the endogenous social effect, and $\beta_3^0$ measures the exogenous social effect. We let $\beta^0 := (\beta^0_1, \beta_2^{0'}, \beta_3^{0'})'$ and denote $\beta = (\beta_1, \beta_2^{'}, \beta_3^{'})'$.
We let $\mathbf{D}_N$ be the $ (N \times N)$ adjacency matrix of the network whose $(i,j)^{th}$ element is $d_{ij}$. We let $d_{ii}=0$ for all $i$, following convention. Let $\mathbf{G}_N$ be the matrix whose $(i,j)^{th}$ element is $g_{ij}$. Recall that $\mathbf{G}_N$ is obtained by row-normalizing $\mathbf{D}_N$. Denote $\mathbf{X}_{1N}=(\mathbf{x}_{11}',\ldots,\mathbf{x}_{1N}')'$, $\mathbf{y}_N=(y_1,\ldots,y_N)'$ and $\bm{\upsilon}_N=(\upsilon_1,\ldots,\upsilon_N)'$. {\color{black} Using this notation, we can express the linear-in-means peer effects model ((ref)) as
Throughout the paper, we assume that $| \beta_1^0 | < 1$. It is known that when $\mathbf{G}_N$ is row normalized (i.e., $\sum_{j \neq i}g_{ij} = 1$) and $ | \beta_1^0 | < 1$, the (equilibrium) solution of the peer effect model uniquely exists (e.g., see Bramoulle2009) as
} In the standard linear-in-means model of peer effects, the main focus has been identification and estimation of peer effects, assuming that the peer group (or the network) is exogenous, that is, $\mathbb{E}[\upsilon_i|\mathbf{X}_{1N},\mathbf{G}_N]=0$. For example, see Manski1993a and Bramoulle2009, Lee2007, and blume2015linear. To identify and estimate the linear-in-means model of peer effects when the peer group is exogenous, it is necessary to take into account the fact that the regressor $\sum_{i=1}^N g_{ij}y_{j}$ is correlated with the error term $\upsilon_i$. For example, if $\upsilon_i\sim\ i.i.d. (0,\sigma^2)$, it is true that
To solve this endogeneity problem different estimators have been proposed in the literature, see for example Kelejian1998, Lee2003 and Lee2007a. One of the widely used estimation methods is the Instrumental Variables (IV) approach. {\color{black} In view of the expression of ((ref)), when $\beta_2^0 \neq 0$, we can use $\mathbf{G}^2_N\mathbf{X}_{1N}$ as the IV of the endogenous regressor $\mathbf{G}_N\mathbf{y}_N$ because $\mathbf{G}^2_N\mathbf{X}_{1N}$ is uncorrelated with $\bm{\upsilon}_N$ while it is correlated with the endogenous regressor $\mathbf{G}_N \mathbf{y}_N$ (see for example Kelejian1998, Lee2003, and Bramoulle2009)\footnote{ {\color{black} If $\beta_2^0 = 0$, $\mathbf{y}_N$ does not depend on $\mathbf{X}_{1N}$ and $\mathbf{G}_N^2\mathbf{X}_{1N}$ is not a relevant instrument for $\mathbf{G}_N \mathbf{y}_N$. }}. }Then, the natural estimator is the Two-Stage Least Squares (2SLS) estimator,
where $\mathbf{W}_N=[\mathbf{G}_N\mathbf{y}_N,\ \mathbf{X}_{1N},\ \mathbf{G}_N\mathbf{X}_{1N}]$ and $\mathbf{Z}_N=[\mathbf{X}_{1N},\ \mathbf{G}_N\mathbf{X}_{1N},\ \mathbf{G}^2_N\mathbf{X}_{1N}]$ is the matrix of instruments. For the IVs $\mathbf{Z}_N$ to be strong, we assume that $\beta_2^0 \neq 0$.
When the network matrix is endogenous, $\mathbb{E}[\mathbf{G}_N\bm{\upsilon}_N]\neq 0$, and the procedure used by Kelejian1998, Lee2003, Bramoulle2009 and others is no longer valid since the IV matrix $\mathbf{Z}_N=[\mathbf{X}_{1N},\ \mathbf{G}_N\mathbf{X}_{1N},\ \mathbf{G}^2_N\mathbf{X}_{1N}]$ is correlated with the error term $\bm{\upsilon}_N$. Specifically, the validity of the 2SLS estimator depends on the orthogonality condition $\mathbb{E}[\bm{\upsilon}_N|\mathbf{Z}_N] = 0$, which is implied if $\mathbb{E}[\bm{\upsilon}_N|\mathbf{X}_{1N},\mathbf{G}_N] = 0$. However, it does not hold if the (row normalized) network $\mathbf{G}_N$ is correlated with $\bm{\upsilon}_N$, which is true if unobserved individual characteristics of $\mathbf{G}_N$ directly influence both link formation and individual outcomes.
In this paper, we consider the case where it may be that $\mathbb{E}[\bm{\upsilon}_N|\mathbf{X}_{1N},\mathbf{G}_N]\neq 0$, so that unobserved characteristics that influence link formation can also have a direct effect on individual outcomes. This is an important consideration in many common applications, like the impact of school friendships on scholarly achievement or substance use. Imagine kids from homes where parents help with homework who only form friendships with kids from similar homes. If this unobserved characteristic of parental behavior is not taken into account, and if this is what really determines grades, this effect might falsely be classified as a peer effect. {\color{black} A more elaborate discussion of our framework and its empirical applications can be found in Section (ref).}
Let $\mathbf{x}_{2i}$ be a vector of observable characteristics of individual $i$, and let $\mathbf{x}_i=\mathbf{x}_{1i}\cup \mathbf{x}_{2i}$. Define $\mathbf{X}_{2N}$ analogously to $\mathbf{X}_{1N}$ and let $\mathbf{X}_N=\mathbf{X}_{1N}\cup\mathbf{X}_{2N}$. We introduce $a_i$, a scalar unobserved characteristic of individual $i$, which is treated as an individual fixed effect, and hence, might be correlated with $\mathbf{x}_i$. We denote the vector of individual unobserved characteristics by $\mathbf{a}_N=(a_1,a_2,\ldots,a_N)'$. Individuals are connected by an undirected network $\mathbf{D}_N$, with the $(i,j)^{th}$ element $d_{ij} = 1$ if $i$ and $j$ are directly connected and $0$ otherwise. We assume the network to be undirected\footnote{{\color{black}Our analysis can be extended to the directed network case, but we do not pursue it in this paper.}}, $d_{ij} = d_{ji}$, and assume $d_{ii}=0$ for all $i$, following the convention. In this case, there are $n=\binom{N}{2}$ dyads. Let $\mathbf{t}_{ij}$ denote an $l_T\times 1$ vector of dyad-specific characteristics of dyad $ij$, and we assume that $\mathbf{t}_{ij}=t(\mathbf{x}_{2i},\mathbf{x}_{2j})$. Agents form links according to
where $\mathbb{I}( \bullet)$ is an indicator function. In this setup, link surplus is transferable across directly linked agents and consists of three components: $\mathbf{t}_{ij} := t(\mathbf{x}_{2i},\mathbf{x}_{2j}) $ is a systematic component that varies with observed dyad attributes and accounts for homophily, $a_i$ and $a_j$ account for unobserved dyad attributes (degree heterogeneity), and $u_{ij}$ is an idiosyncratic shock that is i.i.d. across dyads and independent of $\mathbf{t}_{ij}$ and $a_i$ for all $i,j$. Since links are undirected, the surplus of link $d_{ij}$ must be the same for individual $i$ and $j$. Hence, we assume that the function $t_{ij}$ is symmetric in $i$ and $j$, and the function $g$ is symmetric in $a_i$ and $a_j$.
In the literature, various parametric versions of the network formation in ((ref)) are used, ({\color{black} see for example jackson2005survey, Graham2017)}). An important example of a parametric specification is the one in Graham2017,
For the purpose of the paper, particularly in constructing the estimators that we introduce in Section (ref), we do not need a parametric specification.
Regarding the network formation ((ref)), we impose restrictions (Assumption (ref) (iii) - (vi) in the Appendix) that imply the following two features. The first feature is that the link formation probability of individual $i$ with characteristics $(\mathbf{x}_{2i},a_i)$ is one-to-one with respect to the unobserved characteristic $a_i$, that is, for all $x_{2i}$,
Obviously, this condition is satisfied in the parametric model ((ref)). This monotonic condition justifies the use of the average node degree in implementing the control function as introduced in Section (ref) and will be discussed in Section (ref). The second feature is that the network formed by ((ref)) is dense in the sense that the expected number of connections is proportional to the square of the network size. This is satisfied if the error $u_{ij}$ is drawn randomly from a distribution with full support, while $g( \mathbf{t}_{ij}, a_i, a_j )$ is bounded (see Assumption (ref) (iii),(iv), and (v) in the Appendix). In this case, the probability of any two individuals forming a link is bounded away from zero and strictly less than one. The dense network model is appropriate for scenarios where any two individuals can plausibly form a link. Notice that the dense network assumption and the sharing restriction on the net surplus function $g$ are necessary for implementing the control function in Section (ref) and establishing the asymptotic theory of the control function based estimators in Section (ref). If $a_i$ is observed, we can identify and estimate peer effects without these assumptions (see Section (ref)).
Regarding the network formation model ((ref)), it is important to note that the network formation model ((ref)) rules out interdependent link preferences, and it assumes that links are formed independently conditional on observed individual characteristics and unobserved fixed effects. As discussed in Graham2017, this assumption is appropriate for settings where link formation is driven predominantly by bilateral concerns, such as certain types of friendship networks, trade networks and some models of conflict between nation-states. The model in ((ref)) is not a good choice when important strategic aspects influence link formation, like when the identity of the nodes to which $j$ is linked influences $i$'s return from forming a link with $j$. A discussion of networks with interdependent links can be found in Graham2017 and dePaula2016. Also, when network externalities are present, the additional complication of multiple equilibria has to be considered, see for example Sheng2012 for more details.
In this section we provide an identification argument for the peer effect equation based on a control function when the network is endogenous.
In this subsection we discuss how to control the endogeneity of the peer group defined by the network formed in equation ((ref)). First we introduce a basic assumption that we will maintain throughout the paper.
Assumption (ref)(i) implies that the observables $\mathbf{x}_i$ and the unobservable characteristics $(a_i,\upsilon_i)$ are randomly drawn. This is a standard assumption in the peer effects literature. Assumption (ref)(ii) assumes that the link formation error $u_{ij}$ is orthogonal to all other observables and unobservables in the model. This means that the dyad-specific unobservable shock $u_{ij}$ from the link formation process does not influence outcomes $(y_1,\ldots,y_N)'$. However, we allow for endogeneity of the social interaction group through dependence between the two unobserved components $a_i$ and $\upsilon_i$. This means that the unobserved error $\upsilon_i$ in the outcome equation can be correlated with unobserved individual characteristics $a_i$ that are determinants of link formation. We also allow the observed characteristics $\mathbf{x}_i$ of the outcome equation and the network formation to be correlated with the unobserved components $(\upsilon_i,a_i)$, so that the regressor $\mathbf{x}_{1i}$ can be endogenous in the outcome equation, and the network formation observables $\mathbf{x}_{2i}$ can be arbitrarily correlated with the unobserved individual characteristic $a_i$. In Assumption (ref)(iii), we assume that the dependence between $\mathbf{x}_i$ and $\upsilon_i$ exists only through $a_i$. That is, $a_i$ is the fixed effect of individual $i$ and controls the endogeneity of $\mathbf{x}_i$ with respect to $\upsilon_i$.
Notice that the network $\mathbf{D}_N$ defined in ((ref)) and the (row normalized) network $\mathbf{G}_N$ are measurable functions of $ (\mathbf{x}_{2i},\mathbf{x}_{2,-i},a_i,\mathbf{a}_{-i},\{u_{ij}\}_{i,j=1,\ldots,N}),$ where $\mathbf{x}_{2,-i}=(\mathbf{x}_{2,1},\ldots,\mathbf{x}_{2,i-1},\mathbf{x}_{2,i+1},\ldots,\mathbf{x}_{2,N})$ and $\mathbf{a}_{-i}$ is defined analogously. Under Assumption (ref) we have
where the second equality holds because $(\mathbf{x}_{-i},\mathbf{a}_{-i},\{u_{ij}\}_{i,j=1,\ldots,N})$ and $ (\mathbf{x}_{i},a_i,\upsilon_i)$ are independent under Assumptions (ref) (i) and (ii). This shows $v_i$ and $(\mathbf{x}_{-i}, \mathbf{G}_N(\mathbf{x}_{2,-i},\mathbf{a}_{-i},\{u_{ij}\}_{i,j=1,\ldots,N},\mathbf{x}_{2i},a_i))$ are mean-independent conditioning on $(\mathbf{x}_{i},a_i)$. The last line follows by the fixed effect assumption, Assumption (ref) (iii).
Result ((ref)) shows that conditional on the unobserved heterogeneity $a_i$ in the network formation (and any subcomponents of $\mathbf{x}_i$), the unobserved characteristic $\upsilon_i$ that affects the outcome $y_i$ becomes uncorrelated with the (row normalized) network $\mathbf{G}_N$ (and the observables $\mathbf{X}_N$). This implies that the network endogeneity can be controlled by $a_i$ (or together with any subcomponents of $\mathbf{x}_i$). We summarize the discussion above in the following lemma:
In this section we show how to identify the peer effects in the outcome question when the endogenous network is formed by ((ref)). We provide two identification methods depending on whether we control the network (peer group) endogeneity with $a_i$ or $a_i$ together with $\mathbf{x}_{2i}$, in the case when $\mathbf{x}_{2i}$ and $\mathbf{x}_{1i}$ do not overlap.
First notice that regardless of the possible endogeneity of the (row normalized) network $\mathbf{G}_N$, we need to control for the endogeneity of the term $\sum_{j \neq i} g_{ij}y_j$ that represents the so-called endogenous peer effects. When the peer group $\mathbf{G}_N$ is exogenous and uncorrelated with $\upsilon_N$, $\mathbf{G}^2_N \mathbf{X}_{1N}$ is often used as an IV for the endogenous peer effects term $\mathbf{G}_N \mathbf{y}_N$ (See, for example, Kelejian1998, Lee2003, Bramoulle2009.).
Let $\mathbf{Z}_N=[\mathbf{X}_{1N}, \mathbf{G}_N \mathbf{X}_{1N}, \mathbf{G}^2_N\mathbf{X}_{1N} ]$ be the usual IV matrix used in 2SLS estimation of the peer effects equation. Note that $\mathbf{Z}_N$ is not a valid IV matrix anymore in our framework because the peer group defined by the network $\mathbf{G}_N$ is correlated with $\upsilon_N$ due to potential correlation between the unobserved $\upsilon_i$ and $a_i$. Let $\mathbf{W}_N=[\mathbf{G}_N \mathbf{y}_N, \mathbf{X}_{1N}, \mathbf{G}_N \mathbf{X}_{1N}]$. Further, denote the transpose of the $i$th row of $\mathbf{Z}_N$ and $\mathbf{W}_N$ by $\mathbf{z}_i$ and $\mathbf{w}_i$, respectively.
Suppose that Assumption (ref) holds and so $a_i$ controls the network endogeneity. Then,
where equality $(1)$ holds by Lemma (ref)(a). This shows that the instrumental variables $\mathbf{z}_{i}$ or $\mathbf{z}_{i} -\mathbb{E}[\mathbf{z}_{i}|a_i]$ become orthogonal to $ \upsilon_i-\mathbb{E}[\upsilon_i|a_i],$ the residual of $\upsilon_i$ after projecting out $a_i$.
Furthermore, if $\mathbb{E}\left[ \left( \mathbf{z}_i -\mathbb{E}[\mathbf{z}_i|a_i] \right) \left( \mathbf{w}_i - \mathbb{E}[\mathbf{w}_i|a_i]\right)^{\prime}\right]$ has full rank, then we can identify the peer effect coefficients $\beta^0$ as
where equality $(1)$ follows by the orthogonality result in ((ref)) and equality $(2)$ follows from the full rank condition.
For the full rank condition in Assumption (ref), it is necessary that the IVs $\mathbf{z}_i$ and the regressors $\mathbf{w}_i$ have additional variation after projecting out the control function $a_i$. As shown in the Supplementary Appendix (ref), when $N$ is large, both $\mathbf{z}_i$ and $\mathbf{w}_i$ become close to functions that depend only on $(\mathbf{x}_i,a_i)$. In this case, for the full rank condition to be satisfied, it is necessary that there be additional random components in $\mathbf{x}_{i}$ that are different from $a_i$, so that the limits of $\mathbf{z}_i$ and $\mathbf{w}_i$ are not linearly dependent. As a summary, we have the following first identification theorem.
Theorem (ref) shows that we can identify the parameter $\beta^0$ by controlling the unobserved network heterogeneity $a_i$ in the outcome equation and taking the residuals $y_i - \mathbb{E} (y_i|a_i) - (\mathbf{w}_i-\mathbb{E}(\mathbf{w}_i|a_i))'\beta$ and using the instrumental variables $\mathbf{z}_i - \mathbb{E}[\mathbf{z}_i|a_i]$.
In view of the derivation of the control function in ((ref)) under Assumption (ref), it is possible to use any regressors in $\mathbf{x}_i$ in addition to the unobserved heterogeneity $a_i$. In this section, we discuss identification of the peer effects using $(\mathbf{x}_{2i},a_i)$ as control function. The reason to consider this particular control function is that we can {\color{black} implement it in the absence of a consistent estimator of $a_i$}, which will be discussed in detail in Section (ref).
First, suppose that there is no overlap between the regressors in the outcome equation $\mathbf{x}_{1i}$ and the regressors in the network formation equation $\mathbf{x}_{2i}$ and assume the conditions in Assumption (ref).\footnote{Later in this section, we will discuss a more general case where $\mathbf{x}_{1i}$ and $\mathbf{x}_{2i}$ intersect.}
Then, under Assumption (ref) and by ((ref)), it follows that
where the last line holds by Assumption (ref)(iii). Then, similar to ((ref)), we can show that
Furthermore, suppose that the following full rank assumption is satisfied:
Notice that if $\mathbf{x}_{1i}$ and $\mathbf{x}_{2i}$ are overlapped, then the full rank condition in Assumption (ref) does not hold.
Using similar arguments that lead to Theorem (ref), we can identify the peer effect coefficients $\beta^0$ as
This is summarized in the following theorem.
So far, we have considered the case where the regressors $\mathbf{x}_{i1}$ and $\mathbf{x}_{2i}$ do not intersect. A more general case is when the regressors $\mathbf{x}_{1i}$ consist of two components, where one component is different from the observed control function $\mathbf{x}_{2i}$ and the other is part of $\mathbf{x}_{2i}$. That is, $\mathbf{x}_{1i} = (\mathbf{x}_{11i},\mathbf{x}_{12i})$, where $\mathbf{x}_{11i}$ does not share any elements with $\mathbf{x}_{2i}$ and $\mathbf{x}_{11i}$ is nonempty, and $\mathbf{x}_{12i} \subset \mathbf{x}_{2i}$. Let $\beta^0_2 = (\beta^0_{21},\beta^0_{22}), \beta^0_3 = (\beta^0_{31},\beta^0_{32})$ conformable to the dimensions of $(\mathbf{x}_{11i},\mathbf{x}_{12i})$. Similarly let $\beta_2 = (\beta_{21},\beta_{22}), \beta_3 = (\beta_{31},\beta_{32}).$
In this case, with a properly modified rank condition of $\mathbf{z}_{(2),i}$ and $\mathbf{w}_{(2),i}$ which excludes the variables associated with $\mathbf{x}_{12,i}$ and $\sum_{j=1,\neq i}^N g_{ij} \mathbf{x}_{12,j}$, we can identify the coefficients $\beta^0_{(2)} := (\beta^0_1,\beta^0_{21},\beta^0_{31})$ using the same argument that leads to the identification in ((ref)). However, we cannot identify the coefficients that correspond to the variable $\mathbf{x}_{12,i}$ and $\sum_{j=1,\neq i}^N g_{ij} \mathbf{x}_{12,j}$. The reason is that controlling the network endogeneity with the control variable $(\mathbf{x}_{2i},a_i)$ wipes out the information in $(\mathbf{x}_{12,i},\sum_{j=1,\neq i}^N g_{ij} \mathbf{x}_{12,j})$:
where the second convergence holds because $\sum_{j=1,\neq i}^N g_{ij} \mathbf{x}_{12,j}$ converges to a function that depends only on $(\mathbf{x}_{2i}, a_i)$ (see Section (ref) in the Supplementary Appendix.).
Throughout the rest of the paper, when we consider $(\mathbf{x}_{2i},a_i)$ as control function, we will without loss of generality apply the restriction in Assumption (ref) that $\mathbf{x}_{1i}$ and $\mathbf{x}_{2i}$ do not overlap.
In this section we present two estimation methods. In subsections (ref) and (ref) we discuss estimation using $a_i$ and $(\mathbf{x}_{2i},a_i)$ as control functions, respectively.
The identification scheme of Theorem (ref) identifies the parameter of interest $\beta^0$ with the two step procedure: (i) control $a_i$ in the outcome equation and yield $y_i - \mathbb{E}(y_i|a_i) = (\mathbf{w}_i - \mathbb{E}(\mathbf{w}_i|a_i))'\beta^0 + \upsilon_i - \mathbb{E}(\upsilon_i)$, and then (ii) use $\mathbf{z}_i - \mathbb{E}(\mathbf{z}_i|a_i)$ as IVs for $\mathbf{w}_i - \mathbb{E}(\mathbf{w}_i|a_i)$. If we observe $a_i$ and know the conditional mean functions $\mathbf{h}(a_i)=(h^y(a_i),\mathbf{h}^{w}(a_i),\mathbf{h}^z(a_i)) :=(\mathbb{E}[y_i|a_i], \mathbb{E}[\mathbf{w}_i|a_i],\mathbb{E}[\mathbf{z}_i|a_i])$, then $\beta^0$ can be estimated using 2SLS as
However, since the individual heterogeneity $a_i$ is not observed and the conditional mean functions $\mathbf{h}(a_i) = (\mathbb{E}(y_i|a_i), \mathbb{E}(\mathbf{w}_i|a_i),\mathbb{E}(\mathbf{z}_i|a_i))$ are not known either, the estimator $\widehat{\beta}_{2SLS}^{\text{inf}}$ is not feasible.
A natural implementation of the infeasible estimator $\widehat{\beta}_{2SLS}^{\text{inf}}$ is to replace the conditional mean function $\mathbf{h}(a_i)$ with its estimate. Suppose that $\widehat{a}_i$ is an estimator of $a_i$ and $\widehat{\mathbf{h}}(\widehat{a}_i)$ is a nonparametric estimator of $\mathbf{h}(a_i)$. Then we can implement the infeasible estimator $\widehat{\beta}_{2SLS}^{\text{inf}}$ with
See Section (ref) in the Appendix for more details on the estimator $\widehat{\beta}_{2SLS}$.
{ \bf Estimation of $\mathbf{h}(\cdot)$:} We can estimate $\mathbf{h}(\cdot)$ using various standard nonparametric methods. In this paper we consider a (linear) sieve estimation method.\footnote{In principle we can use other nonparametric estimation methods such as kernel smoothing or local polynomial methods.} Suppose that $h^l(a)$ is the $l^{th}$ element in $\mathbf{h}(a)$ for $l=1,...,L$, where $L$ is the dimension of $(y_i,\mathbf{w}_i',\mathbf{z}_i')'$. The sieve estimation method assumes that each function $h^l(a)$, $l=1,...,L$ is well approximated by a linear combination of base functions $(q_1(a),...,q_{K_N}(a))$:
as the truncation parameter $K_N \rightarrow \infty$. A linear sieve (or series) estimator of a function, for example $\widehat{h}^y(\widehat{a}_i)$, is the OLS projection of $y_i$ on the sieve basis $\mathbf{q}^K(\cdot) = (q_1(\cdot),...,q_K(\cdot))'$ with $\widehat{a}_i$ plugged in, \[ \widehat{h}^y(\widehat{a}_i) := \mathbf{q}^K(\widehat{a}_i)' \left( \sum_{i=1}^N \mathbf{q}^K(\widehat{a}_i)\mathbf{q}^K(\widehat{a}_i)' \right)^{-1} \sum_{i=1}^N \mathbf{q}^K(\widehat{a}_i)y_i. \]
For the regularity conditions of the sieve basis $\mathbf{q}^K(a_i)$, we impose standard conditions such as those proposed by Newey1997 and Li2008a. These assumptions ensure that $\sum_{i=1}^N \mathbf{q}^K(a_i) \mathbf{q}^K(a_i)'$ is asymptotically non-singular and control the rate of approximation of the sieve estimator. These assumptions are formally stated in Assumptions (ref) and (ref) of the Appendix.
Additionally, we require that the sieve basis satisfy a Lipschitz condition, which allows us to control for the error introduced by the estimation of $a_i$ with $\widehat{a}_i$ in the estimation of $\widehat{\beta}_{2SLS}$\footnote{This issue is similar to the two step series estimation problem in Newey2009. Other papers that investigated the problem of nonparametric or semiparametric analysis with generated regressors include Ahn1993, Mammenetal2012, HahnRidder2013, and Escancianoetal2014, for example.} (see Assumptions (ref) and (ref)). As an example, define the polynomial sieve as follows. Let $Pol(K_N)$ denote the space of polynomials on $[-1,1]$ of degree $K_N$, \[ Pol(K_N)=\left\{ \nu_0 + \sum_{k=1}^{K_N}\nu_k a^k, \ a \in [-1,1], \nu_k \in \mathbb{R} \right\}. \] For any $k$ we have \[ \big| a_1^k - a_2^k \big| = k | \tilde{a}^k | | a_1- a_2 | \leq M k | a_1- a_2 | , \] where $\tilde{a} \in [-1,1]$ and $M$ is a finite constant.
In sieve estimations an important issue is choosing the truncation parameter $K_N$. Well-known procedures for selecting $K_N$ are Mallows' $C_P$, generalized cross-validation and leave-one-out cross-validation. For more on these methods see Chapter 15.2 in Li2008a, li1987asymptotic, wahba1985comparison, andrews1991asymptotic and hansen2014nonparametric. However, these methods are mainly applicable when the observations are cross-sectionally independent, which is not true in our case, especially when the network is dense, as we assume. Developing a data-driven choice of $K_N$ is beyond the scope of this paper and we leave it for future work.
{ \bf Estimation of $a_i$:} A desired estimator of $a_i$ should satisfy the following high level condition.
{\color{black} Here $\zeta_a(N)$ is the order of magnitude that measures the Lipschitz smoothness of the sieve basis. The assumption puts restrictions on the uniform bound of the convergence rate of $\widehat{a}_i$, and we need a more accurate estimator of $a_i$ when the average curvature of the sieve basis is larger.}
For the purpose of our paper, any estimation method that yields an estimator $\widehat{a}_i$ satisfying the restriction in Assumption (ref) can be adopted. For example, assuming the parametric specification as in ((ref)),
with regularity conditions of Assumption (ref) in the Appendix, {\color{black} including the error $u_{ij}$ following a logistic distribution,} Graham2017 showed that the joint maximum likelihood estimator that solves
satisfies
with probability $1-O(N^{-2})$. In this case we have $\zeta_a(N) = \sqrt{\frac{N}{\ln N}}$. Notice that the requirement that the network formation in ((ref)) be dense is necessary for $\widehat{a}_i$ to satisfy the desired uniform convergence rate in ((ref)). Examples of other estimation methods include Fernandez-val, jochmans2016modified, dzemski2017empirical, and jochmans2018semiparametric.
As we assume in Section (ref), we consider the case where $\mathbf{x}_{1i}$ and $\mathbf{x}_{2i}$ do not overlap. When $a_i$ is observed and the conditional expectations $\mathbf{h}_{*}(\mathbf{x}_{2i},a_i) = (h_{*}^y(\mathbf{x}_{2i},a_i),\mathbf{h}_{*}^w(\mathbf{x}_{2i},a_i),\mathbf{h}_{*}^z(\mathbf{x}_{2i},a_i)):= (\mathbb{E}(y_i|\mathbf{x}_{2i},a_i),\mathbb{E}(\mathbf{w}_i|\mathbf{x}_{2i},a_i),\mathbb{E}(\mathbf{z}_i|\mathbf{x}_{2i},a_i))$ are known, we can estimate $\beta^0$ by the 2SLS similar to $\widehat{\beta}^{\inf}_{2SLS}$ in ((ref)),
When $a_i$ is unknown and $\mathbf{x}_{2i}$ is also used in the control function, under the monotonicity condition of the link formation as in ((ref)), we can implement the infeasible estimator using the average node degree without estimating $a_i$. To be more specific, first we denote
Under the monotonicity condition in ((ref)), $(\mathbf{x}_{2i},a_i)$ and $(\mathbf{x}_{2i}, \text{deg}_i)$ are one-to-one. This implies that for any $b_i \in \{ y_i,\mathbf{w}_{i},\mathbf{z}_{i} \}$, \[ h_*^b(\mathbf{x}_{2i},a_i) = \mathbb{E}(b_i| \mathbf{x}_{2i},a_i) = \mathbb{E}(b_i| \mathbf{x}_{2i},{\rm deg}_i) =: h_{**}^b(\mathbf{x}_{2i},{\rm deg}_i). \]
Notice that the natural estimator of ${\rm deg}_i$ is the node degree of $i$, the number of connections with node (individual) $i$ in the network scaled by the network size: \[ \widehat{\text{deg}}_i := \frac{1}{N-1} \sum_{j=1, \neq i}^N d_{ij}. \] Recall that the link $d_{ij}$ is formed by \[ d_{ij}=\mathbb{I}( g( t(\mathbf{x}_{2i},\mathbf{x}_{2j}), a_i, a_j ) - u_{ij} \geq 0). \] Also recall that the unobserved link-specific error terms $u_{ij}$ are assumed to be independent of all the other variables and randomly drawn. Let $\Phi(\cdot)$ be the cdf of $u_{ij}$. Also let $\pi(\mathbf{x}_2,a)$ be the joint density function of $(\mathbf{x}_{2i},a_i)$. Then, for each $(\mathbf{x}_{2i},a_i)$, by the WLLN conditioning on $(\mathbf{x}_{2i},a_i)$, we have
as the network size $N$ grows to infinity. Here the limit of the average network $\rm{deg}_i > 0$ follows since we assume the network is dense.
This shows that $\widehat{\text{deg}}_i$ can be used as an estimator of $\text{deg}_i$. In fact, we can show that under the regularity conditions in Assumption (ref) in the Appendix, $\sup_{i} \mathbb{E} [( \sqrt{N} (\widehat{\text{deg}}_i - \text{deg}_i ))^{2B} ] < \infty$ for any finite integer $B \geq 2$, from which we can deduce that
where
This corresponds to the regularity condition in Assumption (ref).
Suppose that $\mathbf{r}^{K}(\mathbf{x}_{2i},\text{deg}_i) =(r_1(\mathbf{x}_{2i},\text{deg}_i),\ldots,r_{K}(\mathbf{x}_{2i},\text{deg}_i))'$ is a sieve basis of the unknown function $\mathbf{h}_*(\mathbf{x}_{2i},a_i)$. For each $b_i \in \{ y_i,\mathbf{w}_{i},\mathbf{z}_{i} \}$, a sieve estimator of $h_{**}^b(\mathbf{x}_{2i},\text{deg}_i) = \mathbb{E}(b_i| \mathbf{x}_{2i},a_i)$ is the OLS projection of $b_i$ on $\mathbf{r}^{K}(\mathbf{x}_{2i},\widehat{\text{deg}}_i)$. For example,
Then, we have
For more details see Section (ref) in the Appendix.
The two different estimators $\widehat{\beta}_{2SLS}$ and $\bar{\beta}_{2SLS}$ are implemented using different control functions, and these two approaches have their own pros and cons. For $\widehat{\beta}_{2SLS}$, a good estimator of $a_i$ is required, which imposes restrictions on the network formation model ((ref)) in the form of ((ref)). Compared to this, the estimator $\bar{\beta}_{2SLS}$ that uses $(\mathbf{x}_{2i}, {\rm deg}_i)$ as control functions does not require a restriction like ((ref)). It requires only the monotonicity of the net surplus function as in ((ref)) of Section (ref). However, $\bar{\beta}_{2SLS}$ has disadvantages: because it uses $x_{2i}$ as a part of the control function, as discussed in Section (ref), this approach cannot identify and estimate the coefficients of the regressor $\mathbf{x}_{2i}$ if $\mathbf{x}_{2i}$ is a relevant regressor of the outcome. Later in Section (ref), where we present the Monte Carlo simulations, we compare the finite sample properties of $\widehat{\beta}_{2SLS}$ and $\bar{\beta}_{2SLS}$ in both dense and sparse network setups.
In this section we present the asymptotic distributions of the two 2SLS estimators $\widehat{\beta}_{2SLS}$ and $\bar{\beta}_{2SLS}$, and show how to estimate standard errors. We also discuss key technical issues in deriving the limits. All details of the technical derivations and proofs can be found in the Appendix.
Recall the definitions $h^{y}(a_i):= \mathbb{E}[y_i|a_i], \quad h^{\upsilon}(a_i):= \mathbb{E}[\upsilon_i|a_i], \quad \mathbf{h}^{\mathbf{w}}(a_i) := \mathbb{E} (\mathbf{w}_i|a_i), \quad \mathbf{h}^{\mathbf{z}}(a_i) := \mathbb{E} (\mathbf{z}_i|a_i).$ Define $\eta^{y}_i: = y_i - h^{y}(a_i), \quad \eta^{\upsilon}_i: = \upsilon_i - h^{\upsilon}(a_i), \quad \eta^{\mathbf{w}}_i = \mathbf{w}_i - \mathbf{h}^{\mathbf{w}}(a_i), \quad \eta_i^{\mathbf{z}} = \mathbf{z}_i - \mathbf{h}^{\mathbf{z}}(a_i).$ Let $\bm{\eta}_N^{\upsilon} = (\eta^{\upsilon}_1,...,\eta^{\upsilon}_N)'$ and $\mathbf{H}^{\upsilon}_N(\mathbf{a}_N) = (h^{\upsilon}(a_1),...,h^{\upsilon}(a_N))'$. Let $\widehat{h}^{\upsilon}(a_i)$, $\widehat{\mathbf{h}}^{\mathbf{w}}(a_i)$, and $\widehat{\mathbf{h}}^{\mathbf{z}}(a_i)$ denote the sieve estimators of $h^{\upsilon}(a_i)$, $h^{\mathbf{w}}(a_i)$ and $h^{\mathbf{z}}(a_i)$, respectively.
In the Appendix, we derive the asymptotic distribution of $\widehat{\beta}_{2SLS}$ in three steps. First, we show that the sampling error caused by the use of $\hat{a}_i$ instead of $a_i$ is asymptotically negligible (see Lemma (ref) of the Supplementary Appendix (ref).). Next, we control the error introduced by the non-parametric estimation of $h^{l}(a_i)$, where $l \in \{\upsilon,\mathbf{w},\mathbf{z}\}$. In Lemma (ref) of the Supplementary Appendix (ref) we show that under the regularity conditions, the estimation error in $\widehat{h}^l(a_i)$ vanishes at a suitable rate. Combining these two, we deduce \[ \sqrt{N} (\widehat{\beta}_{2SLS} - \widehat{\beta}^{\inf}_{2SLS}) = o_p(1). \] The last step is to derive the limiting distribution of the infeasible estimator $\sqrt{N} ( \widehat{\beta}^{\inf}_{2SLS} - \beta^0)$. In the Supplementary Appendix (ref) we show the following:
where the closed forms of the limits $\mathbf{S}^{\mathbf{w}\mathbf{z}}$ and $\mathbf{S}^{\mathbf{z}\mathbf{z}}$ are found in Lemma (ref) and $\mathbf{S}^{\mathbf{z}\mathbf{z}\sigma}$ in Lemma (ref) of Supplementary Appendix.
Notice that the derivation of the limiting distribution in ((ref)) allows $\eta^{\upsilon}_i = \upsilon_i - \mathbb{E}(\upsilon_i|a_i)$ to be conditionally heteroskedastic, and so $\sigma^2(\mathbf{x}_i,a_i) := \mathbb{E}[(\upsilon_i - \mathbb{E}[ \upsilon_i|a_i])^2| \mathbf{x}_i,a_i]$ is allowed to depend on $(\mathbf{x}_i,a_i)$.
Combining all the limit results leads to the following theorem.
The theorem requires several regularity conditions which are presented in Appendix (ref). In addition to conditions of random sampling of $(y_i,\mathbf{x}_i,a_i)$ in Assumption (ref) and the full rank condition in Assumption (ref), we assume conditions that ensure $a_i$ can be consistently estimated, and that the error between $\mathbf{h}(a_i)$ and $\widehat{\mathbf{h}}(\widehat{a}_i)$ converges to zero at a suitable rate (Assumptions (ref), (ref) and (ref)). We also impose restrictions on the outcome model ((ref)) and the network formation model ((ref)) (Assumption (ref)). We assume $|\beta_1^0|$ is bounded below $1$ so that the spillover effect has a unique solution, and $\| \beta_2^0 \|$ is bounded above $0$ so that the IVs are strong. We also assume the observables $(y_i,\mathbf{x}_i)$ and $\mathbf{t}_{ij}$ are bounded, and $a_i$ has a compact support in $[-1,1]$. This boundedness condition is required as a technical regularity condition that simplifies the proofs of the limits in ((ref)), ((ref)), and ((ref)), which involves some uniformity in the limit.
The asymptotic variance can be consistently estimated by
where
and $\widehat{\eta}^{\upsilon}_i=y_i-\widehat{h}^y(\widehat{a}_i)-(\mathbf{w}_i-\widehat{\mathbf{h}}^{\mathbf{w}}(\widehat{a}_i))'\widehat{\beta}_{2SLS}.$
The process is analogous to the one presented in the previous section. Again, let $b_i^l$ be the $l^{th}$ element in $(y_i,\mathbf{w}_i',\mathbf{z}_i')'$. Recall the definition that
Further, let $\eta^l_{*i}=b^l_i-h^l_{*}(\mathbf{x}_{2i},a_i) = b^l - h^l_{**}(\mathbf{x}_{2i},\text{deg}_i)$, and let $\widehat{h}^l_{**}(\mathbf{x}_{2i},\text{deg}_i)$ denote a sieve estimator of $h^l_{**}(\mathbf{x}_{2i},\text{deg}_i)$.
As in the previous section, we derive the asymptotic distribution of $\bar{\beta}_{2SLS}$ in three steps. First, we show that the error that stems from the use of the estimate $\widehat{\text{deg}_i}$ for $\text{deg}_i$, $\widehat{h}^l_{**}(\mathbf{x}_{2i},\widehat{\text{deg}}_i) - \widehat{h}^l_{**}(\mathbf{x}_{2i},\text{deg}_i)$, is asymptotically negligible. In the second step, we control the error introduced by the non-parametric estimation of $h^l_{**}(\mathbf{x}_{2i},\text{deg}_i)$, $\widehat{h}^l_{**}(\mathbf{x}_{2i},\text{deg}_i)-h^l_{**}(\mathbf{x}_{2i},\text{deg}_i)$. This implies \[ \sqrt{N} (\bar{\beta}_{2SLS} - \bar{\beta}^{\inf}_{2SLS}) = o_p(1). \] The last step is to derive the limiting distribution of the infeasible estimator $\sqrt{N} ( \bar{\beta}^{\inf}_{2SLS} - \beta^0)$ by showing
Combining all the limit results we have the following theorem.
The asymptotic result in Theorem (ref) requires the following regularity conditions which are formally presented in the Appendix. First, Assumption (ref) assumes that the regressors in the outcome equation, $\mathbf{x}_{1i}$ and the observables in the network formation $\mathbf{x}_{2i}$ do not overlap. Assumption (ref) is a full rank condition for $\bar{\beta}_{2SLS}$. Assumptions (ref) and (ref) regard the sieve used in constructing the estimator $\bar{\beta}_{2SLS}$. Comparing with the assumptions assumed in Theorem (ref), Theorem (ref) does not require the high level condition of Assumption (ref) because we do not use an estimator of $a_i$. Instead it requires an additional restriction that the net surplus function in the link formation be strictly monotonic in $a_i$ conditional on $(\mathbf{x}_{2i},\mathbf{x}_{2j},a_j)$, which implies the required monotonicity condition in ((ref)).
Like in the case of $\widehat{\beta}_{2SLS}$, we allow $\eta^{\upsilon}_{*i} = \upsilon_i - \mathbb{E}(\upsilon_i|\mathbf{x}_{2i},a_i)$ to be conditionally heteroskedastic, and $\sigma^2_{*}(\mathbf{x}_i,a_i) := \mathbb{E}[(\upsilon_i - \mathbb{E}[ \upsilon_i|\mathbf{x}_{2i},a_i])^2| \mathbf{x}_i,a_i]$ is allowed to depend on $(\mathbf{x}_i,a_i)$.
The asymptotic variance can be consistently estimated by
where
and $\widehat{\eta}_{**i}^{\upsilon}=y_i-\widehat{h}_{**}^y(\mathbf{x}_{2i},\widehat{\text{deg}}_i)-(\mathbf{w}_i-\widehat{\mathbf{h}}_{**}^{\mathbf{w}}(\mathbf{x}_{2i},\widehat{\text{deg}}_i))'\bar{\beta}_{2SLS}.$
We consider both dense and sparse network Monte Carlo designs. In the dense network case links are formed according to\footnote{This follows the approach of Graham2017.} \[ d_{ij} = \mathbb{I}\left\{ x_{2i}x_{2j}\lambda_d + a_i + a_j -u_{ij} \geq 0 \right\}, \] where $x_{2i}\in\{-1,1\}$, $\lambda_d=1$ and $u_{ij}$ follows a logistic distribution. This link rule implies that agents have a strong taste for homophilic matching since $x_{2i}x_{2j}\lambda_d=1$ when $x_{2i}=x_{2j}$ and $x_{2i}x_{2j}\lambda_d=-1$ when $x_{2i}\neq x_{2j}$.
In the sparse network case links are formed according to \[ d_{ij} = \mathbb{I}\left\{(|x_{2i}-x_{2j}|+3)\lambda_s + a_i + a_j -u_{ij} \geq 0 \right\}, \] with $\lambda_s=-1$. This rule also implies homophily on observable characteristics. Individual-level degree heterogeneity is generated according to \[ a_i=\varphi(\alpha_L\mathbb{I}\left\{ x_{2i}=-1 \right\} +\alpha_H\mathbb{I} \left\{ x_{2i}=1 \right\} + \xi_i), \] with $\alpha_L \leq \alpha_H$ and $\xi_i$ a centered Beta random variable $ \xi_i|x_{2i}\sim \left\{Beta(\mu_0,\mu_1)-\frac{\mu_0}{\mu_0+\mu_1}\right\}$ so that $a_i\in \left[\alpha_L-\frac{\mu_0}{\mu_0+\mu_1},\alpha_H+\frac{\mu_1}{\mu_0+\mu_1}\right]$. We choose values of the network formation parameters so that $a_i \in[-1,1]$. In the main text we present results based on the following parameter values. In the dense network case we set $\mu_0=1/4$, $\mu_1=3/4$, $\alpha_L=\alpha_H=-3/4$, which yields an average node degree $=23$ when $N=100$. The sparse network formation design is generated by setting $\mu_0=1$, $\mu_1=1$, $\alpha_L=\alpha_H=-1/4$, which gives an average degree $=1.78$ when $N=100$.\footnote{Results for 14 other network formation designs can be found in Section (ref) of the online appendix. Most results are similar to the ones presented in the main text.}
Individual outcomes are generated according to \[ y_i=\beta_1\sum_{j=1 \atop j\neq i}^N g_{ij}y_j+\beta_2 x_{1i}+\beta_3 \sum_{j=1 \atop j\neq i}^N g_{ij}x_{1j}+h(a_i)+\varepsilon_i. \] In the simulations, we set $\beta_1=0.8$, $\beta_2=\beta_3=5$, $x_{1i}=3q_1+\cos(q_2)/0.8+\epsilon_i$, where $q_1,q_2\sim\mathcal{N}(x_{2i},1)$, and $\varepsilon_i,\epsilon_i\sim\mathcal{N}(0,1)$. For $h(a_i)$ we use the following functional forms: $h(a_i)=\exp (3 a_i)$, $h(a_i)=\cos(3 a_i)$, $h(a_i)=\sin(3 a_i)$. A plot of $h(a_i)$ for these functional forms is presented in Figure (ref). We can see that the exponential function yields a strongly increasing impact on the individual outcome, and with the cosine functions the returns are increasing up to a certain point and then decreasing; however the sine function gives a more irregular pattern.
We estimate the outcome equation coefficients $(\beta_1,\beta_2, \beta_3)$ using the standard 2SLS estimator for peer effects and the Hermite polynomial sieve as well as a polynomial sieve. For the dense network case, we estimate $a_i$ using $\widehat{a}_i$ and implement the following control functions: using a control function linear in $\widehat{a}_i$, $\widehat{h}(\widehat{a}_i)$, $\widehat{h}(a_i)$, $\widehat{h}(\widehat{\rm deg}_i,x_{2i})$\footnote{ {Note that since $x_{2i}$ is discrete with a finite support, $\{ x_1,...,x_M \}$,} we have $ r(x_{2i},{\rm deg_i}) = \sum_{m=1}^M r(x_m,{\rm deg_i}) \mathbb{I}\{ x_{2i} = x_m \}. $ We can then approximate $ r(x_{2i}, {\rm deg_i}) \simeq \sum_{k=1}^{K_N} \left\{ \sum_{m=1}^M \alpha_{m,k} q_k^d(\rm deg_i) \mathbb{I}\{ x_{2i} = x_m \} \right\}. $}, and $h(a_i)$. For the sparse network case the estimator of $a_i$ is not reliable\footnote{To estimate $a_i$, we use the JMLE proposed in Graham2017. As Graham2017 states, in sparse designs the JMLE rarely even exists, rendering it unusable in practice when the network is too sparse. See Graham2017 for more details.} and we implement the following control functions: linear in $a_i$, $\widehat{h}(a_i)$, $\widehat{h}(\widehat{\rm deg}_i,x_{2i})$ and $h(a_i)$. In both the dense and sparse setup we also implement a benchmark model with no control for the endogeneity of the network.
In the paper, due to space limitations, we present Monte Carlo results obtained using the Hermite polynomial sieve with $K_N=4$. Specifically, Tables (ref) and (ref) include results for the dense and sparse network specifications, respectively. Results for the other orders of $K_N$ are not notably different; in the Online Supplement we provide results for fourteen other network formation designs, for $K_N=4,8$ and for the Hermite polynomial and polynomial sieve functions.
We also perform conventional leave-one-out cross validation to find data-dependent $K_N$ (chosen as the $K_N$ that minimizes the Root Mean Square Error (RMSE) of the prediction based on the leave-one-out estimator, see for example andrews1991asymptotic, hansen2014nonparametric). We report the statistics on the cross-validation in Table (ref). The differences in RMSE are very small between the different values of $K_N$.\\ Analyzing the Monte Carlo results for the dense network specification in Table (ref), we can see that, as expected from our asymptotic theories, the control functions $\widehat{h}(\widehat{a}_i)$ and $\widehat{h}(\widehat{deg}_i,x_{2i})$ perform better than the estimator with a linear control function, as well as the estimator that does not control for the endogeneity of the network in terms of mean bias. This difference is more pronounced in the case when $h(a_i)$ is the sine or cosine function. Both the control for degree approach and the control function that uses $\widehat{h}(\widehat{a}_i)$ yield a low bias and have the correct size on all coefficients in all cases. In the simulations we also implemented the control function $\widehat{h}(a_i)$, that is, using the true $a_i$ instead of $\widehat{a}_i$. These results are very similar to the ones obtained using $\widehat{h}(\widehat{a}_i)$, which is in line with the estimator $\widehat{a}_i$ having a very low bias, as detailed in the table footnotes. This suggest that the approach of using $\widehat{h}(\widehat{a}_i)$ as a control function works very well when a highly precise estimator of $a_i$ is available (for example when the network size $N$ is large.).
Looking at Table (ref) and the results for the sparse design, we can see that the control for degree approach performs very well across all functional forms of $h(a_i)$. In the sparse setup, the bias of all estimates, including those that do not control for the endogeneity of the network, is small. However, the size of the no control and linear control estimates is not correct. If a precise estimator of $a_i$ is available, the control function $\widehat{h}(a_i)$ also performs well with low bias and correct size in all cases.
Table (ref) shows that the performance of the estimators does not differ notably for different values of $K_N$. As for the choice of $K_N$ we present in the tables, we have run simulations for a range of values of $K_N$ and the results did not differ significantly. As deriving a theory for a data driven choice of $K_N$ is beyond the scope of this paper, for applied researchers we suggest estimating the model over a range of $K_N$ and seeing whether the results vary significantly. As shown in our Monte Carlo simulations, the control function approach yields results robust to the choice of $K_N$ for different non-linear functions.
In this paper we show that, whenever the network is likely endogenous, it is important to control for this endogeneity when estimating peer effects. Failing to control for the endogeneity of the connections matrix in general leads to biased estimates of peer effects. We show that under specific assumptions, we can use the control function approach to deal with the endogeneity problem. We assume that unobserved individual characteristics directly affect link formation and individual outcomes. We leave the functional form through which unobserved individual characteristics enter the outcome equation unspecified and estimate it using a non-parametric approach. The estimators we propose are easy to use in applied work, and Monte Carlo results show that they perform well compared to a linear control function estimator. Erroneously assuming that unobserved characteristics enter the outcome equation in a linear fashion can lead to a serious bias in the estimated parameters.
{
}