$\acute{\tilde{Y}}_{iOLS}(\hat{\beta}_{t},\hat{\Lambda}_{t})=(I_{n}-P_{0})\tilde{Y}_{iOLS}(\hat{\beta}_{t},\hat{\Lambda}_{t})$
and the remaining regressors variable $\acute{X_{1}}=(I_{n}-P_{0})X_{1}$;

  4. 4.

Compute the OLS estimate
$\hat{\beta}_{t+1}=[\acute{X_{1}}^{\prime}\acute{X_{1}}]^{-1}\acute{X_{1}}^{\prime}\acute{\tilde{Y}}(\hat{\beta}_{t})$,
and update $t$ to $t+1$;

  5. 5.

Recover the fixed-effects into the aggregate term
$\hat{\Lambda}_{t}=(\tilde{Y}(\hat{\beta}_{t})-X_{1}^{\prime}\hat{\beta}_{t+1})-(\acute{\tilde{Y}}(\hat{\beta}_{t})-\acute{X_{1}}^{\prime}\hat{\beta}_{t+1})$

  6. 6.

Iterate steps 2 to 5 until $\hat{\beta}_{t}$ converges.

Note that all matrix inversions in this algorithm can be done only once. The
presence of fixed-effects has hence almost no effect on the computation speed
of the iterative estimator. Remark further that this approach relates to the
Poisson estimator with high-dimensional fixed-effects. Indeed, Correia,
Guimarães and Zylkin (2019) transform the PPML estimator into an iteratively
reweighted least squares problem, then make use of the Frisch-Waugh-Lovell
theorem to speed up computations like above. Their approach bears some
similarities with our approach for iOLSAP (additive poisson), except that ours
involves naturally less matrix inversions. The i2SLS estimator extends
similarly to many fixed-effects.

##### Incidental parameter problem.

In non-linear panel data models, individual fixed-effects are not always
consistent when the cross-sectional dimension $n$ increases to infinity while
the time dimension $T$ remains fixed. This issue is known as the incidental
parameters problem (IPP). It is a well-known issue with maximum likelihood
estimators, and some solutions have been recently developed (Fernández-Val and
Weidner, 2016; Weidner and Zylkin, 2021). Although we do not study this
problem from a theoretical viewpoint, our replications of the simulations in
Weidner and Zylkin (2021) reveal that the iOLS estimators exhibit an IPP bias
very similar to PPML. Adapting the solutions proposed in the literature should
hence be possible. We leave this topic for future research.

###  B.6 The log of a ratio

Researchers are sometimes willing to estimate equations of the form

$\log(Y_{i1}/Y_{i2})=X_{i}^{\prime}\beta+\varepsilon_{i},$ (112)

where $Y_{i1}$ and $Y_{i2}$ are two outcomes of interest. It may happen that
both outcomes can take zero values, hence not only the log is undefined but
also the ratio. The “popular fix” estimator in this case consists in
transforming the outcomes and focus on the following model

$\log((Y_{i1}+\Delta)/(Y_{i2}+\Delta))=X_{i}^{\prime}\beta+\omega_{i},$ (113)

for some $\Delta>0$.Needless to explain why this simple fix is not
satisfactory. Instead, let us consider an alternative solution where the two
following equations are estimated jointly

$\displaystyle\log(Y_{i1}+\Delta)$
$\displaystyle=X_{i}^{\prime}\beta_{1}+\varepsilon_{1i}$ (114)
$\displaystyle\log(Y_{i2}+\Delta)$
$\displaystyle=X_{i}^{\prime}\beta_{2}+\varepsilon_{2i},$

by rewriting the problem as a seemingly unrelated regression problem. Here, we
propose to use the popular fix as a starting point, but other methods like
iOLS will apply without difficulty. The seemingly unrelated regression model
can be easily implemented as a single iOLS regression. The parameter $\beta$
of interest corresponds to $\beta_{1}-\beta_{2}$ and inference can be
conducted using the delta-method. The advantage of this approach is that one
can separately check which model is best to address the log of zero in each
equation.

###  B.7 An alternative iOLS transformation

An alternative iOLS transformation would consist in letting $\delta$ vary
across observations. For example, let $\delta_{i}=\delta(1-\xi_{i})$ where
$\xi_{i}$ takes a zero value when $Y_{i}=0$ and is equal to 1 otherwise.
Therefore, the iOLS transform becomes

$\log(Y_{i}+(1-\xi_{i})\delta\exp(X_{i}^{\prime}\beta))=X^{\prime}\beta+\log((1-\xi_{i})\delta+U_{i}).$
(115)

Let us recall that $U_{i}=\exp(\varepsilon_{i})\xi_{i}$, thus the error term
is $\log((1-\xi_{i})\delta+\exp(\varepsilon_{i})\xi_{i})$. We now develop its
conditional mean into

$E(\log((1-\xi_{i})\delta+\exp(\varepsilon_{i})\xi_{i})|X)=E(\varepsilon_{i}|\xi_{i}=1,X)P(X)+\log(\delta)(1-P(X))$
(116)

On the other hand, the exogeneity condition imposed in the log-linear model is
about

$E(\varepsilon_{i}|X)=E(\varepsilon_{i}|\xi_{i}=1,X)P(X)+E(\varepsilon_{i}|\xi_{i}=0,X)(1-P(X)).$
(117)

Therefore, imposing the restriction
$E(\log((1-\xi_{i})\delta+\exp(\varepsilon_{i})\xi_{i})|X)=0$ under the
assumption that $E(\varepsilon_{i}|X)=0$ (log-linear) is equivalent to
assuming that

$E(\varepsilon_{i}|\xi_{i}=0,X)=\log(\delta),$ (118)

where $\delta$ can be chosen using the testing procedures presented in the
paper.

More generally, $E(\varepsilon_{i}|X)=0$ implies that

$E(\varepsilon_{i}|\xi_{i}=1,X)=-E(\varepsilon_{i}|\xi_{i}=0,X)(1-P(X))P(X)^{-1}.$
(119)

We can hence evaluate any assumption about $E(\varepsilon_{i}|\xi_{i}=0,X)$ by
considering a function $\delta(\cdot)>0$ and test whether the following
condition holds

$E(\varepsilon_{i}|\xi_{i}=1,X)=-\log(\delta(X))(1-P(X))P(X)^{-1}.$ (120)

This approach can be helpful although the choice of the candidate functions
for $\delta(\cdot)$ to be considered is beyond the scope of this paper. Remark
that Heckman’s selection model corresponds to specifying

$\delta(X)=\exp(-\lambda\frac{\phi(-X_{i}^{\prime}\gamma)}{1-\Phi(X_{i}^{\prime}\gamma)}).$
(121)

###  B.8 Testing with endogenous regressors

In this section, we explain how our tests adapt to endogenous regressors.

##### Testing the Poisson condition.

For Poisson models, we have

$E[U_{i}|Z_{i}]=E[U_{i}|Z_{i},U_{i}>0]Pr(U_{i}>0|Z_{i})=E(U_{i}),$ (122)

since $E[U_{i}|Z_{i},U_{i}=0]=0$. Following the same step as with exogenous
regressors, the error term $U_{i}$ under the null is such that

$U_{i}=\lambda E[U]Pr(U_{i}>0|Z_{i})^{-1}+\nu_{i}$ (123)

for $U_{i}>0$ with $\lambda=1$ and $E[\nu_{i}|U_{i}>0,Z_{i}]=0$. There are
hence two differences: 1. one needs to estimate $P(U>0|Z)$ instead of
$P(U>0|X)$, and 2. an IV estimator, like i2SLS, must be used to obtain
$\hat{U}$.

##### Testing the i2SLS restriction.

For iOLSδ, we have $E[\log(\delta+U_{i})|Z_{i}]=c$. The null hypothesis is now

$H_{0}:E[\log(\delta+U_{i})|Z_{i},U_{i}>0]-\log(\delta)=\frac{c-\log(\delta)}{Pr(U_{i}>0|Z_{i})},$
(124)

hence the differences are the same than for Poisson models.

##### Testing other restrictions.

Testing for other restrictions introduces some new steps. Developing the
associated exogeneity condition yields

$E[\omega_{i}|Z_{i},U_{i}>0]P(Z_{i})+E[\omega_{i}|Z_{i},U_{i}=0](1-P(Z_{i}))=0$
(125)

which can be rearranged into

$E[\omega_{i}|Z_{i},U_{i}>0]=-E[\omega_{i}|Z_{i},U_{i}=0](1-P(Z_{i}))P(Z_{i})^{-1}.$
(126)

For the popular fix estimator, substituting the expression of $\omega_{i}$ on
the RHS gives

$E[\omega_{i}|Z_{i},U_{i}>0]=-(\log(\Delta)-E(X^{\prime}\beta|Z,U>0))(1-P(Z_{i}))P(Z_{i})^{-1},$
(127)

where the new term $E(X^{\prime}\beta|Z,U>0)$ can be obtained from the first-
stage estimates of the 2SLS procedure neglecting the zero values. For the IHS
estimator, we have the similar form

$E[\omega_{i}|Z_{i},U_{i}>0]=E(X^{\prime}\beta|Z,U>0)(1-P(Z_{i}))P(Z_{i})^{-1}.$
(128)

##  Appendix C Additional simulations results

###  C.1 iOLS DGP

We also considered a simulation where the “true DGP” is such that
$E[\log(\delta+U_{i})-c)|X]=0$, as required by iOLSδ. This DGP is useful to
illustrate iOLSδ under ideal conditions. We fix $\delta=0.25$ and $c=-0.7447$
so that $E(U)\approx 1$. We assume that $(X_{1i},X_{2i})^{\prime}$ is like in
DGP 1. We further assume that $\log(\delta+\exp(\varepsilon))$ follows a
truncated Gaussian distribution with mean
$\frac{c-\log(\delta))}{P(X)}+\log(\delta)$, and with $\log(\delta)$ as lower
bound to guarantee $\exp(\varepsilon)>0$. The log-scale error $\varepsilon$ is
hence approximately Gaussian. The variance parameter is fixed at 0.5. Note
that the error $U_{i}$ is also heteroskedastic in this DGP.

Table C.1 reports the results for DGP 3 (iOLS). All estimators but iOLSδ are
biased. The automatic selection of $\delta$ does not seem to introduce much
noise in the smaller sample. Interestingly, the precision of PPML and iOLSAP
estimates does not seem to improve with the sample size in this DGP.

Table C.1:  Simulations: DGP 3 (iOLS)

|  | $n=1000$ |  |  | $n=10,000$ |
---|---|---|---|---|---|---
Estim. | $\beta_{0}$ | $\beta_{1}$ | $\beta_{2}$ | $\beta_{0}$ | $\beta_{1}$ | $\beta_{2}$
$iOLS_{\delta}$ (best) | $0.71$ | $1.02$ | $0.98$ | $0.84$ | $1.00$ | $1.00$
| $(0.33)$ | $(0.13)$ | $(0.13)$ | $(0.28)$ | $(0.05)$ | $(0.05)$
$iOLS_{\delta}$ (auto) | $0.78$ | $1.00$ | $1.00$ | $0.88$ | $1.00$ | $1.00$
| $(0.62)$ | $(0.13)$ | $(0.12)$ | $(0.45)$ | $(0.04)$ | $(0.04)$
$iOLS_{MP}$ | $0.44$ | $0.67$ | $1.33$ | $0.47$ | $0.63$ | $1.37$
| $(0.17)$ | $(0.11)$ | $(0.11)$ | $(0.10)$ | $(0.08)$ | $(0.08)$
$iOLS_{AP}$ | $0.29$ | $0.43$ | $1.46$ | $0.08$ | $0.24$ | $1.66$
| $(1.10)$ | $(0.45)$ | $(0.47)$ | $(1.40)$ | $(0.46)$ | $(0.51)$
$PPML$ | $0.30$ | $0.43$ | $1.46$ | $0.08$ | $0.24$ | $1.66$
| $(1.12)$ | $(0.46)$ | $(0.48)$ | $(1.38)$ | $(0.46)$ | $(0.51)$
$OLS$ | $1.09$ | $0.51$ | $1.49$ | $1.09$ | $0.51$ | $1.49$
| $(0.08)$ | $(0.05)$ | $(0.05)$ | $(0.02)$ | $(0.01)$ | $(0.01)$
$IHS$ | $0.66$ | $0.53$ | $0.26$ | $0.66$ | $0.53$ | $0.26$
| $(0.10)$ | $(0.07)$ | $(0.07)$ | $(0.03)$ | $(0.02)$ | $(0.02)$
$PF$ | $0.49$ | $0.46$ | $0.26$ | $0.49$ | $0.46$ | $0.26$
| $(0.08)$ | $(0.06)$ | $(0.06)$ | $(0.03)$ | $(0.02)$ | $(0.02)$

  * •

Notes: This table shows the mean and standard errors (in parentheses) of
parameter estimates across 10,000 simulations based on DGP3.

###  C.2 KNN-based specification tests

As explained in Section 4.1, the properties of the tests depend on the
specified conditional probability function. To address this concern, one can
opt for using a non-parametric estimator of the conditional probability. For
comparison, we report the same statistics when using kNN instead of the
parametric logit in Table C.2. We select a 100 neighbors for the estimation.
Our choice is based on the observation that the larger the number of
neighbors, the smoother the resulting probability estimates because they are
allowed to take values in a larger set, here 101 possible values. Fewer number
of neighbors might possibly lead to better classification predictions but
introduces more noise in the estimation of $\lambda$. The results show that
the non-parametric estimator introduces some distortions in the tests’ sizes
but exhibit reasonably good power.

Table C.2:  Simulations: Specification testing (kNN)

DGP |  |  iOLSδ | MP | PPML | RESET
---|---|---|---|---|---
| $\delta$ | 0.01 | 0.15 | 1.9 | 24 | 86 |  |  |
1 | $\lambda$ | $1.02$ | $1.02$ | $1.01$ | $1.01$ | $1.01$ | $1.00$ | $1.01$ |
(se) | $(0.01)$ | $(0.01)$ | $(0.01)$ | $(0.01)$ | $(0.01)$ | $(0.01)$ | $(0.08)$ |
Rej% | $92.8$ | $78.1$ | $25.4$ | $6.9$ | $9.7$ | $3.2$ | $6.2$ | $7.3$
2 | $\lambda$ | $1.01$ | $1.00$ | $0.98$ | $0.93$ | $0.91$ | $0.87$ | $Inf$ |
(se) | $(0.01)$ | $(0.01)$ | $(0.01)$ | $(0.02)$ | $(0.03)$ | $(0.06)$ | $(Inf)$ |
Rej% | $47.8$ | $13.1$ | $44.9$ | $93.3$ | $91.8$ | $97.3$ | $20.1$ | $6.3$
3 | $\lambda$ | $0.99$ | $0.97$ | $0.92$ | $0.82$ | $0.76$ | $0.66$ | $Inf$ |
(se) | $(0.04)$ | $(0.06)$ | $(0.09)$ | $(0.10)$ | $(0.11)$ | $(0.15)$ | $(Inf)$ |
Rej | $11.0$ | $28.2$ | $92.3$ | $98.6$ | $98.0$ | $96.7$ | $0.8$ | $5.9$

  * •

Notes: This table shows the relative rejection frequency of each null
hypothesis for 1,000 simulations for $n=10,000$ and DGP 1 to 3.

###  C.3 Misspecification of the conditional probability

For DGPs 1 to 3, we have also used alternative specifications of the
conditional probability function to investigate how it affects the automatic
selection of $\delta$. Figure C.1 shows that the value of the selected
$\delta$ is not much affected by the specified function.

(a) DGP1 (Poisson)

(b) DGP2 (iOLS)

(c) DGP3 (Log-linear)

Figure C.1: iOLS’s $\hat{\lambda}$ as functions of $\delta$ with different
specification of the conditional probability function (n=10,000)

###  C.4 Endogenous regressors and two-way fixed-effects

We consider two additional DGPs specified as follows:

  * •

DGP 1 with endogenous regressors (IV): $E[U_{i}|X_{i}]\neq 1$ but
$E[U_{i}|Z_{i}]=1$. Let us assume that
$Pr(\xi_{i}=0|Z_{i})=P(Z_{i})=\frac{1}{1+\exp(\gamma_{0}+\gamma_{1}Z_{1i}+\gamma_{2}Z_{2i})}$,
with the same parameters as in the main text. The instrumental variables
$Z_{1i}$ and $Z_{2i}$ are iid normal with mean 1 and variance
$\sigma_{Z1}^{2}=\sigma_{Z2}^{2}=1$. We further assume that $\varepsilon_{i}$
is Gaussian with mean $-\log(P(Z_{i}))-1/2$ and variance 1 so that
$\exp(\varepsilon_{i})$ is log-normal with conditional mean $1/P(Z_{i})$.
Finally the endogenous regressors are such that
$X_{ik}=0.8Z_{ik}+0.2\varepsilon_{i}^{2}$, for $k=1,2$.

  * •

DGP 1 with endogenous regressors and fixed-effects (IV-FE): the DGP is the
same as above except for individual and time fixed-effects, denoted
$\alpha_{i}$ and $\rho_{t}$, and assumed to be iid uniformly distributed in
$[-0.5,0.5]$. There are $n/T$ individuals with $T=100$ periods, for a total of
$1+T+n/T$ regressors.

Table C.3 reports the results for DGP IV, where the regressors are endogenous
and requires the use of instrumental variables to achieve identification under
the assumption that $E[U_{i}|Z_{i}]=1$. The i2SLS estimators perform similarly
to the non-linear IV estimator (Mullahy, 1997) ($IV-MP$) in terms of bias and
variance. However, the control function approach of Wooldridge (1997) ($CF-
MP$) delivers biased estimates because the endogenous regressors are non-
linear functions of the error term in this design. The estimators denoted $IV-
AP$ and $CF-AP$ correspond to the same estimators based on the additive
Poisson model. They fail to deliver consistent estimates. $2SPF$ denotes the
two-stage popular fix, where $\Delta=1$ is added to $Y$ before taking the log
function in the second stage.

i2SLS estimators bring important computational gains in settings with fixed-
effects. We only simulate this model 10 times to evaluate the computational
complexity of each estimator. For $n=1,000$, hence 111 regressors,
$i2SLS_{MP}$ is computed in about 0.08 second whereas $IV-MP$ takes 45
seconds. For $n=10,000$, hence 201 regressors, $i2SLS_{MP}$ is computed in
about 0.8 second whereas $IV-MP$ takes 500 seconds.

Table C.3:  Simulations: DGP 1 with endogenous regressors (DGP IV)

|  | $n=1000$ |  |  | $n=10,000$ |
---|---|---|---|---|---|---
Estim. | $\beta_{0}$ | $\beta_{1}$ | $\beta_{2}$ | $\beta_{0}$ | $\beta_{1}$ | $\beta_{2}$
$i2SLS_{\delta}$ (best) | $1.07$ | $0.99$ | $0.97$ | $1.00$ | $1.01$ | $0.99$
| $(1.28)$ | $(0.38)$ | $(0.16)$ | $(0.08)$ | $(0.04)$ | $(0.03)$
$i2SLS_{\delta}$ (auto) | $1.07$ | $1.15$ | $0.84$ | $1.00$ | $1.03$ | $0.96$
| $(0.37)$ | $(0.25)$ | $(0.25)$ | $(0.08)$ | $(0.05)$ | $(0.04)$
$i2SLS_{MP}$ | $1.08$ | $0.98$ | $0.98$ | $1.00$ | $1.00$ | $1.00$
| $(1.25)$ | $(0.38)$ | $(0.16)$ | $(0.09)$ | $(0.05)$ | $(0.03)$
$IV-MP$ | $1.04$ | $0.99$ | $0.98$ | $1.00$ | $1.00$ | $1.00$
| $(0.36)$ | $(0.17)$ | $(0.10)$ | $(0.09)$ | $(0.05)$ | $(0.03)$
$IV-AP$ | $-0.88$ | $1.37$ | $1.19$ | $0.18$ | $1.24$ | $1.10$
| $(14.09)$ | $(3.51)$ | $(2.11)$ | $(12.10)$ | $(2.14)$ | $(1.00)$
$CF-MP$ | $0.74$ | $1.08$ | $0.92$ | $0.74$ | $1.07$ | $0.93$
| $(0.13)$ | $(0.08)$ | $(0.07)$ | $(0.04)$ | $(0.03)$ | $(0.02)$
$CF-AP$ | $-Inf$ | $Inf$ | $Inf$ | $1.32$ | $1.18$ | $0.78$
| $(Inf)$ | $(Inf)$ | $(Inf)$ | $(2.07)$ | $(0.73)$ | $(0.62)$
$2SPF$ | $0.62$ | $0.69$ | $0.07$ | $0.62$ | $0.69$ | $0.07$
| $(0.14)$ | $(0.09)$ | $(0.07)$ | $(0.04)$ | $(0.03)$ | $(0.02)$

  * •

Notes: This table shows the bootstrapped parameter estimates and standard
errors calculated on data simulated according to DGP IV, as described in
Section 5 and Appendix C. The column “Estim.” reports the different estimated
models. The estimates are reported based on a sample of size $n=1000$ or of
$n=10,000$. Standard Errors are presented in between parentheses and are
calculated using pairs bootstrap based on 10,000 simulations.

##  Appendix D Data Appendix

###  D.1 American Economic Review (2016-2020)

Table D.1: Solutions to the Log of Zero in the AER (2016-2020)

Log of Zero |  log($\Delta$ \+ Yi) |  PPML |  Drop |  IHS
---|---|---|---|---
48 |  23 (48%) |  17 (35%) |  15 (31%) |  7 (15%)

  * •

Notes: This table reports the number of articles published in the American
Economic Review from 2016 to 2020 where the issue of the log of zero was
encountered. “Log of Zero” is the number of publications where at least one
regression had to address this issue. “log($\Delta$ \+ Yi)” refers to the
common fix of adding some discretionary constant to the dependent variable
before taking the logarithmic transformation. “PPML” refers to Pseudo-Poisson
Maximum Likelihood or Negative Binomial regression. “Drop” refers to cases
where the problematic observations are discarded. “IHS” refers to the Inverse
Hyperbolic Sine Transformation of the dependent variable. Some articles used
several solutions, as robustness checks, which explains why the sum of
solutions is different larger than 48.

Table D.2: American Economic Review Cases per Year  Year |  Emp. Pub. |  log(Yi) |  log($\Delta$+Yi) |  PPML |  Drop |  IHS
---|---|---|---|---|---|---
2016 |  69 |  27 |  2 |  4 |  7 |  1
2017 |  71 |  28 |  5 |  2 |  4 |  1
2018 |  69 |  32 |  4 |  4 |  2 |  1
2019 |  79 |  27 |  6 |  6 |  2 |  3
2020 |  53 |  19 |  6 |  1 |  0 |  1

  * •

Notes: This table displays the frequency of solutions observed in American
Economic Review. The sample extends over the period Jan. 2016 to Oct. 2020.
Emp. Pub. is the number of empirical papers (includes “data” section). The
column $\log(Y_{i})$ counts cases where the dependent variable was in
logarithmic form or in which a fix (such as $\log(\Delta+Y_{i})$, PPML, Drop,
or IHS) is used. It excludes cases where the author openly states that a
logarithmic specification was preferred but rejected due to the existence of
non-positive observations. $\log(\Delta+Y_{i})$ is the popular fix. PPML
refers to Poisson and Negative Binomial regression. Drop refers to cases where
the author dropped the problematic observations. IHS is the Inverse Hyperbolic
Transformation.

###  D.2 ResearchGate

Figure D.1: Proposed solutions by category on ResearchGate (November 2018)

###  D.3 Wooclap Survey

Figure D.2: Wooclap Survey

Description. The survey was implemented during 3 seminars (CREST, HEC
Montréal, and University of Montréal) in 2021, before the speaker presented
the different approaches. The attendees could provide multiple answers to the
questions displayed in Figure D.2 and were invited to indicate if they were a
student. Results are presented in Table D.3.

Table D.3: Wooclap Survey Results

| Frequency
---|---
Popular fix | 42,8 %
Poisson | 17,8 %
Other transformation | 17,8 %
Mixture | 35,7 %
Drop zeros | 17,8 %
Levels instead of logs | 17,8 %
Another method | 3,5 %
None satisfactory | 25 %
Not applicable | 3,5 %
PhD Student | 21,4 %
Nb. Respondents | 28

  * •

Notes: This table displays relative frequency of answers to the Wooclap
Survey. Intrepretation: 42.8% of respondents would use the popular fix (but
not necessarily exclusively).

###  D.4 Michalopoulos and Papaioannou (2013)

Table D.4: Probability Model of Non-Zero Light Density at Night

Logit Probability Model
---
Dependent var. | $\mathbbm{1}(\text{Light Density at Night}>0)$
Jurisdictional Hierarchy | 0.137 | 0.196
| (0.137) | (0.167)
Ln(Population Density) |  0.690∗∗∗ |  1.078∗∗∗
| (0.112) | (0.215)
Distance to Capital City | -0.429 | -0.745
| (0.562) | (1.045)
Distance to Sea Coast | -0.191 | 1.197
| (0.566) | (1.097)
Distance to Border | -1.513 | -0.577
| (3.112) | (2.625)
Ln(1+Water Area) | -0.0918 | 0.346
| (0.314) | (0.450)
Ln(Land Area) |  1.058∗∗∗ |  1.274∗∗∗
| (0.130) | (0.180)
Mean Elevation | -0.715 | -0.549
| (0.465) | (0.769)
Soil Suitability | -0.415 | -0.0535
| (0.952) | (0.806)
Ecological Suitability |  -1.784∗ | -2.138
| (0.806) | (1.098)
Petroleum Field Dummy | 0.933 | 1.093
| (0.688) | (0.711)
Diamond Mine Dummy | 0.394 | 0.443
| (0.504) | (0.573)
Country Fixed Effect | No | Yes
$N$ | 682 | 621
∗ $p<0.05$, ∗∗ $p<0.01$, ∗∗∗ $p<0.001$

  * •

Notes: This table provides the logit estimates and standard errors (s.e),
based on the research of Michalopoulos and Papaioannou (2013). Each
observation of the data is at the ethnicity by country level. The dependent
variable is a dummy equal to one if the recorded light density at night is
non-zero. Standard errors are provided in between parenthesis are clustered at
the ethnicity by country level.

###  D.5 Santos Silva and Tenreyro (2006)

Table D.5: Replication of Table (3)

| PF | PPML | PPML | PPML |  iOLSMP |  iOLSδ
---|---|---|---|---|---|---
|  |  | (No Zero) | (Top 10%) |  |
Constant | -39.909 | -32.326 | -31.530 | -28.096 | -36.523 | -38.707
| (1.225) | (2.135) | (2.234) | (2.739) | (2.997) | (2.431)
Ln(GDP) - (X) | 1.128 | 0.732 | 0.721 | 0.650 | 0.978 | 1.113
| (0.012) | (0.026) | (0.026) | (0.029) | (0.035) | (0.028)
Ln(GDP) - (I) | 0.866 | 0.741 | 0.732 | 0.672 | 0.859 | 0.951
| (0.012) | (0.028) | (0.029) | (0.034) | (0.027) | (0.023)
Ln(GDP/Capita) - (X) | 0.277 | 0.157 | 0.154 | 0.125 | 0.297 | 0.285
| (0.018) | (0.052) | (0.051) | (0.053) | (0.033) | (0.031)
Ln(GDP/Capita) - (I) | 0.217 | 0.135 | 0.133 | 0.118 | 0.254 | 0.234
| (0.018) | (0.045) | (0.045) | (0.045) | (0.040) | (0.034)
Ln(Distance) | -1.151 | -0.784 | -0.776 | -0.686 | -1.477 | -1.589
| (0.036) | (0.062) | (0.062) | (0.063) | (0.080) | (0.072)
Contiguity | -0.241 | 0.193 | 0.202 | 0.298 | 0.204 | 0.001
| (0.209) | (0.105) | (0.107) | (0.113) | (0.352) | (0.335)
Common Language | 0.742 | 0.746 | 0.751 | 0.702 | 0.692 | 0.877
| (0.065) | (0.138) | (0.138) | (0.136) | (0.132) | (0.126)
Colonial Ties | 0.392 | 0.025 | 0.020 | 0.001 | 0.390 | 0.448
| (0.070) | (0.159) | (0.160) | (0.161) | (0.138) | (0.129)
Landlocked - (X) | 0.106 | -0.863 | -0.872 | -0.886 | -0.589 | -0.494
| (0.056) | (0.153) | (0.153) | (0.150) | (0.130) | (0.114)
Landlocked - (I) | -0.278 | -0.696 | -0.703 | -0.673 | -0.881 | -0.973
| (0.052) | (0.131) | (0.131) | (0.135) | (0.126) | (0.110)
Remoteness - (X) | 0.526 | 0.660 | 0.647 | 0.625 | 1.020 | 0.888
| (0.089) | (0.142) | (0.143) | (0.154) | (0.178) | (0.169)
Remoteness - (I) | -0.109 | 0.561 | 0.549 | 0.563 | 0.113 | 0.036
| (0.090) | (0.125) | (0.126) | (0.133) | (0.218) | (0.162)
Free-trade agreement Dummy | 1.289 | 0.181 | 0.179 | 0.245 | 1.549 | 1.425
| (0.121) | (0.096) | (0.097) | (0.098) | (0.502) | (0.345)
Openness | 0.739 | -0.107 | -0.139 | -0.169 | -0.414 | -0.369
| (0.050) | (0.128) | (0.129) | (0.136) | (0.102) | (0.095)
$N$ | 18,360 | 18,360 | 9,613 | 1,836 | 18,360 | 18,360

  * •

Notes: This table replicates Table 3 in Santos Silva and Tenreyro (2006).
Standard errors based on 300 pairs bootstrap are in parenthesis, and $\lambda$
for various models of trade. iOLSδ and iOLSMP are defined in Section 3. PF is
the baseline relying on the popular fix ($\Delta=1$). PPML (No Zero) denotes
the estimates from dropping zero flows. PPML (Top 10%) denotes the estimates
from keeping the top 25% of trade flows. (X) refers to exporter
characteristics and (I) to importer characteristics.


# Graph Mixer Networks

Ahmet Sarıgün
Middle East Technical University
<EMAIL_ADDRESS>
###### Abstract

In recent years, the attention mechanism has demonstrated superior performance
in various tasks, leading to the emergence of GAT and Graph Transformer models
that utilize this mechanism to extract relational information from graph-
structured data. However, the high computational cost associated with the
Transformer block, as seen in Vision Transformers, has motivated the
development of alternative architectures such as MLP-Mixers, which have been
shown to improve performance in image tasks while reducing the computational
cost. Despite the effectiveness of Transformers in graph-based tasks, their
computational efficiency remains a concern. The logic behind MLP-Mixers, which
addresses this issue in image tasks, has the potential to be applied to graph-
structured data as well. In this paper, we propose the Graph Mixer Network
(GMN), also referred to as Graph Nasreddin Nets (GNasNets), a framework that
incorporates the principles of MLP-Mixers for graph-structured data. Using a
PNA model with multiple aggregators as the foundation, our proposed GMN has
demonstrated improved performance compared to Graph Transformers. The source
code is available publicly at https://github.com/asarigun/GraphMixerNetworks.

##  1 Introduction

Graph Neural Networks (GNNs) are a powerful tool for working with graph-
structured data, which is data that is made up of entities and their
relationships. GNNs have been used to solve a wide range of problems, such as
node classification, link prediction, graph generation, and many others. They
have attracted great interest in recent years due to their performance and the
ability to extract complex information.

Graph Convolutional Networks [1] (GCN) is a type of graph neural network (GNN)
that uses graph convolutional layers to process data represented as graphs.
GCNs can be used for various tasks such as node classification [2], graph
classification [3], and link prediction. In each graph convolutional layer,
the node features are updated by aggregating the features of their neighboring
nodes. This is done through a convolution operation, where a linear
combination of the neighboring node features is applied to each node, followed
by a non-linear activation function. Graph Isomorphism Network (GIN) [4] is
another type of GNN. GIN are able to distinguish non-isomorphic graphs. GIN
consists of multiple layers of neural networks, where each layer aggregates
the features of the neighboring nodes using a sum pooling operation, followed
by a multi-layer perceptron (MLP). GIN can be used for various tasks such as
node classification, graph classification, and link prediction. Graph
Attention Networks (GAT) [5] is a GNN that uses attention mechanisms to assign
different importance to different neighboring nodes when aggregating their
features. Each node in GAT has a self-attention mechanism that allows it to
weigh its own features and the features of its neighboring nodes in a
learnable way. GAT can be used for various tasks such as node classification,
graph classification, and link prediction. Message Passing Neural Networks
(MPNN) [6] is a class of GNNs that generalize the idea of message passing
between nodes in a graph. In MPNNs, messages are passed between nodes in the
graph, and the node updates its state based on the messages received from its
neighbors. MPNNs can be used for various tasks such as node classification,
graph classification, and link prediction.

GCN, GIN, GAT, and MPNN are types of GNNs, each with their own characteristics
and capabilities. GCN uses graph convolutional layers, GIN uses a sum pooling
operation, GAT uses attention mechanisms to assign importance to different
neighboring nodes, and MPNNs generalize the idea of message passing between
nodes in a graph. They can all be used for various tasks such as node
classification, graph classification, and link prediction.

Transformer networks have been particularly successful in natural language
processing tasks such as machine translation, text summarization, and language
understanding.

The key innovation of transformer networks is the use of self-attention
mechanisms. In a transformer network, each element in a sequence (such as a
word in a sentence) is processed by attending to all the other elements in the
sequence. This allows the network to weigh the importance of different
elements in the sequence when making a prediction.

Vision transformers [7] are a type of transformer network designed for
computer vision tasks, adapting the transformer architecture from natural
language processing to handle image data. The input image is divided into non-
overlapping patches and processed individually using self-attention
mechanisms. The key advantage of vision transformers is their ability to
handle images of arbitrary size, unlike traditional convolutional neural
networks [8] which require a fixed image size.

Although transformers have shown superior performance on the image task, they
have a computational quadratic complexity. MLP-Mixers [9, 10] have been used
for the first time in the image task, and it has been shown to perform well
without the quadratic complexity like in the attention mechanism and using
only MLP.

In this work, we propose Graph Mixer Networks (or Graph Nasreddin Networks)
which uses MLP-Mixer on Graphs as a novel Graph Neural Network. They are shown
to have less computational complexity than transformers and comparable
performance to baseline models.

The source code is available publicly at
https://github.com/asarigun/GraphMixerNetworks.

##  2 Graph Mixer Networks

The proposed Graph Mixer Networks (or Graph Nasreddin Networks) method
leverages the increased expressivity from the multi-aggregators models such as
PNA [11] and MLP-Mixer [9, 10]. GMN architecture is given in Figure 1.

Figure 1: Mixer Layer of Graph Mixer Network where $T$ denotes transpose
operator.

###  2.1 Motivation

Graph Transformers use the attention mechanism. Within this mechanism, in the
self-attention mechanism, the dot product between Query and Key makes the
computational complexity by $O(n^{2})$. For this reason, MLP-Mixers, which are
applied in image tasks, are seen to be more computationally efficient and
perform better than Vision Transformers [9, 10]. That is why MLP-Mixer has
considered solving computational efficiency in graphs in this study.

###  2.2 Multi-aggregators

Multiple aggregations are generalizations of sum aggregation within single
aggregators and have been shown theoretically and empirically to be better at
discriminating graphs. Therefore, in this study, we use multiple aggregators
such as min, max, and mean aggregators used in the PNA model simultaneously.

#### Degree Scalers.

In GMN, degree scalers amplify or decrease signals based on a node’s degree,
allowing for more flexibility in the model. The formula for this is $S$
(scaling factor) = $\alpha$ (amplification factor) x $d$ (node degree) /
$delta$ (average degree in the training set). In our research, we tested
amplification factors of $-1$, $0$, and $1$, which respectively result in
attenuation, no change, and signal amplification based on the node’s degree.

$S(d,\alpha)=\left(\frac{\log(d+1)}{\delta}\right)^{\alpha},d>0,-1<\alpha<1$
(1)

#### Combining Aggregators.

Like in the PNA, we incorporate various aggregators and degree scalers using
the equation $\otimes$ (Tensor product) and $\oplus$ (general aggregation
function) in the GMN framework to enhance the network’s flexibility.

$\oplus=\begin{pmatrix}I\\\ S(D,\alpha=1)\\\
S(D,\alpha=-1)\end{pmatrix}\otimes\begin{pmatrix}Max\\\ Min\\\
Mean\end{pmatrix}$ (2)

###  2.3 MLP-Mixer

In this work, we used the Mixer block recommended in MLP-Mixer, which is used
for the image task, as the Mixer layer in GMN. Input,$x$ first goes through
the layer norm operation, then transposes this matrix by the $T$ operator. The
transposed $x$ passes through the Linear layer as $x^{T}$ and is again
transposed and transformed into $x$. It is then merged with the input $x$ with
the residual connection. Then it passes through Layer Norm, linear layer, and
residual connection again. MLP-Mixer block can be summarized as the following
equation 3:

$Mix=MLP_{2}(LayerNorm((MLP_{1}((LayerNorm(x))^{T}))^{T}+x))+x$ (3)

###  2.4 Graph Mixer Network

As in PNA, after $h_{i}^{k}$ and $h_{j}^{k}$, which are the features of the
neighboring node, are concatenated and pass through the linear layer
containing learnable parameters, they pass through multiple aggregation and
scalars. We take the embedding features obtained as a result of this process
here as $x$. This $x$ enters the Mixer block, just like the MLP-Mixers
performed in the image task. Figure 1 shows the components of the Mixer
Blocks, and Equation 3 shows the Mix operator. as a result, the update
function $h_{i}^{k+1}$ is obtained as in Equation 4.

$h_{i}^{k+1}=Mix(h_{i}^{k},\oplus(h_{i}^{k},h_{j}^{k},e_{ij}^{k}))$ (4)

##  3 Experiments

We evaluated the performance of GMN models on ZINC [3] dataset. The
performance results of GMN were compared with the Message Passing Neural
Networks (MPNN) [6] and Graph Transformer.

###  3.1 Dataset

Our method was trained using the ZINC dataset, which is a dataset for
predicting the solubility of chemical compounds through graph regression. The
compounds in the dataset are represented as graphs, with atoms as nodes and
bonds between atoms as edges. The ZINC dataset includes 12,000 molecules, with
atom numbers ranging from 9 to 37. The performance of the method was evaluated
using the mean absolute error (MAE) metric.

###  3.2 Results

We trained models using the multiple aggregator(s) such as mean, max, min.
Finally, we compared our results with the well-known baseline methods in the
literature. The results are given in Table LABEL:table:benchmarking. Our
results have shown improved performance over attention based mechanism
methods.

Table 1: Benchmarking GMN on ZINC dataset. Models |  ZINC(MAE)
---|---
GCN [1] | 0.367
GAT [5] | 0.384
MPNN [6] | 0.288
Graph Transformer [12] | 0.226
GMN (ours) | 0.212
GraphGPS [13] | 0.070

##  4 Discussion and Conclusion

The experiment shows that it is possible to create a powerful transformer-
style graph regressor without using attention layers. Additionally, the MLP-
Mixer model has a significant advantage over the Graph Transformer in terms of
complexity [9, 10], as it is linear in relation to the sequence length instead
of quadratic. This is achieved through the use of an intermediate projection
dimension within the feed-forward layer applied to aggregated learnable
embeddings.

The performance of the GMN is lower than GraphGPS because the positional
encoding [14] and the number of parameters are quite low. In addition to its
lower performance, the main disadvantage of the MLP-Mixer model is that it can
only operate on sequences of a fixed length (as a result of the feed-forward
layer applied to aggregated learnable embeddings). While this is not a problem
in the image domain, it can be a limitation for graph neural networks because
graphs do not have a fixed data structure.

This study shows that MLP-Mixers are effective for graph regression. Future
research should focus on understanding the specific roles of other parts of
the MLP-Mixer, such as interpretability or initialization scheme.
Additionally, the report hopes to inspire further investigation into the
underlying reasons for the performance of current models.

## Acknowledgements

The author express his gratitude to Ahmet S. Rifaioğlu, Gökhan Özsarı and
Mehmet Volkan Atalay not just only for writing and encouragin him to write
this paper but also giving or the valuable insights and discussions.

## References

  * Kipf and Welling [2016] Thomas N Kipf and Max Welling.  Semi-supervised classification with graph convolutional networks.  _arXiv preprint arXiv:1609.02907_ , 2016.
  * Sen et al. [2008] Prithviraj Sen, Galileo Namata, Mustafa Bilgic, Lise Getoor, Brian Galligher, and Tina Eliassi-Rad.  Collective classification in network data.  _AI magazine_ , 29(3):93–93, 2008.
  * Irwin et al. [2012] John J Irwin, Teague Sterling, Michael M Mysinger, Erin S Bolstad, and Ryan G Coleman.  Zinc: a free tool to discover chemistry for biology.  _Journal of chemical information and modeling_ , 52(7):1757–1768, 2012.
  * Xu et al. [2018] Keyulu Xu, Weihua Hu, Jure Leskovec, and Stefanie Jegelka.  How powerful are graph neural networks?  _arXiv preprint arXiv:1810.00826_ , 2018.
  * Veličković et al. [2017] Petar Veličković, Guillem Cucurull, Arantxa Casanova, Adriana Romero, Pietro Lio, and Yoshua Bengio.  Graph attention networks.  _arXiv preprint arXiv:1710.10903_ , 2017.
  * Gilmer et al. [2017] Justin Gilmer, Samuel S Schoenholz, Patrick F Riley, Oriol Vinyals, and George E Dahl.  Neural message passing for quantum chemistry.  In _International conference on machine learning_ , pages 1263–1272. PMLR, 2017.
  * Dosovitskiy et al. [2020] Alexey Dosovitskiy, Lucas Beyer, Alexander Kolesnikov, Dirk Weissenborn, Xiaohua Zhai, Thomas Unterthiner, Mostafa Dehghani, Matthias Minderer, Georg Heigold, Sylvain Gelly, et al.  An image is worth 16x16 words: Transformers for image recognition at scale.  _arXiv preprint arXiv:2010.11929_ , 2020.
  * He et al. [2016] Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun.  Deep residual learning for image recognition.  In _Proceedings of the IEEE conference on computer vision and pattern recognition_ , pages 770–778, 2016.
  * Tolstikhin et al. [2021] Ilya O Tolstikhin, Neil Houlsby, Alexander Kolesnikov, Lucas Beyer, Xiaohua Zhai, Thomas Unterthiner, Jessica Yung, Andreas Steiner, Daniel Keysers, Jakob Uszkoreit, et al.  Mlp-mixer: An all-mlp architecture for vision.  _Advances in Neural Information Processing Systems_ , 34:24261–24272, 2021.
  * Melas-Kyriazi [2021] Luke Melas-Kyriazi.  Do you even need attention? a stack of feed-forward layers does surprisingly well on imagenet.  _arXiv preprint arXiv:2105.02723_ , 2021.
  * Corso et al. [2020] Gabriele Corso, Luca Cavalleri, Dominique Beaini, Pietro Liò, and Petar Veličković.  Principal neighbourhood aggregation for graph nets.  _arXiv preprint arXiv:2004.05718_ , 2020.
  * Dwivedi and Bresson [2020] Vijay Prakash Dwivedi and Xavier Bresson.  A generalization of transformer networks to graphs.  _arXiv preprint arXiv:2012.09699_ , 2020.
  * Rampášek et al. [2022] Ladislav Rampášek, Mikhail Galkin, Vijay Prakash Dwivedi, Anh Tuan Luu, Guy Wolf, and Dominique Beaini.  Recipe for a general, powerful, scalable graph transformer.  _arXiv preprint arXiv:2205.12454_ , 2022.
  * Dwivedi et al. [2021] Vijay Prakash Dwivedi, Anh Tuan Luu, Thomas Laurent, Yoshua Bengio, and Xavier Bresson.  Graph neural networks with learnable structural and positional representations.  _arXiv preprint arXiv:2110.07875_ , 2021.