[a,b]Constantia Alexandrou

# Nucleon transverse quark spin densities

Simone Bacchio  Martha Constantinou  Petros Dimopoulos  Jacob Finkenrath
Roberto Frezzotti  Kyriakos Hadjiyiannakou  Karl Jansen  Bartosz Kostrzewa
Giannis Koutsou  Gregoris Spanoudes  Carsten Urbach

###### Abstract

We present a calculation of the Mellin moments of the nucleon transverse quark
spin densities extracted from the unpolarized and transversity generalized
form factors. We use three $N_{F}=2+1+1$ ensembles of twisted mass fermions
with quark masses tuned to their physical values and lattice spacings $a\sim
0.08$ fm, $a\sim 0.07$ fm and $a\sim 0.06$ fm and extrapolate the form factors
to the continuum limit. Besides isovector densities we also include results
for the tensor charge for each quark flavor using the ensemble with $a\sim
0.08$ fm for which we include the disconnected contributions.

##  1 Introduction

Understanding the 3D-structure of the nucleon from its fundamental
constituents, the quarks and the gluons is a major goal of nuclear physics and
a key aim of current and future major experiments. Significant progress has
been made in revealing the longitudinal spin structure of the proton. However,
the transverse spin structure is much less known and will be the focus of the
experimental programs of Jefferson Lab and the future Electron Ion
Collider(EIC) being built at Brookhaven [1]. Although recently techniques to
compute the generalized parton distributions (GPDs) using lattice QCD are
being developed, computation of moments of GPDs is at a more mature stage. In
particular, lower moments are readily accessible from lattice QCD by computing
matrix elements of local operators. Although the computation of the low Mellin
moments has a long history [2], it is only recently that results directly at
the physical point are emerging (i.e. using ensembles generated with pion mass
at the physical point defined as $m_{\pi}=135\pm 10$ MeV). To determine the
transverse spin densities, one needs to compute the twist-two chiral-even
unpolarized and chiral-odd transversity GPDs.

##  2 Nucleon matrix elements

The transverse spin density is defined as [3]

$\displaystyle\rho(x,{\bf b}_{\perp},{\bf s}_{\perp},{\bf
S}_{\perp})=\frac{1}{2}\bigg{[}H(x,b^{2}_{\perp})+\frac{{\bf
b}_{\perp}^{j}\epsilon^{ji}}{m_{N}}\left({\bf
S}_{\perp}^{i}E^{\prime}(x,b^{2}_{\perp})+{\bf
s}_{\perp}^{i}\bar{E}^{\prime}_{T}(x,b^{2}_{\perp})\right)$ (1)
$\displaystyle+$ $\displaystyle{\bf s}_{\perp}^{i}{\bf
S}_{\perp}^{i}\left(H_{T}(x,b^{2}_{\perp})-\frac{\Delta_{b_{\perp}}\widetilde{H}_{T}(x,b^{2}_{\perp})}{4m_{N}^{2}}\right)+{\bf
s}_{\perp}^{i}(2{\bf b}_{\perp}^{i}{\bf
b}_{\perp}^{j}-\delta^{ij}b^{2}_{\perp}){\bf
S}_{\perp}^{j}\frac{\widetilde{H}^{\prime\prime}_{T}(x,b^{2}_{\perp})}{m_{N}^{2}},$

where the generalized parton distributions $H,E,H_{T},E_{T}$ and
$\widetilde{H}_{T}$ in Eq. (1) are given in impact parameter space for zero
skewness. They are obtained by a Fourier transformation,
${\bf\Delta}_{\perp}\rightarrow{\bf b}_{\perp}$, ${\bf\Delta}_{\perp}$ is the
transverse momentum transfer and $-Q^{2}\equiv\Delta^{2}$. $m_{N}$ is the
nucleon mass, $x$ is the longitudinal momentum fraction, ${\bf s}_{\perp}$ is
the transverse quark spin, ${\bf b}_{\perp}$ the transverse vector from the
center of momentum of the nucleon, and ${\bf S}_{\perp}$ the transverse spin
of the nucleon, $\epsilon_{ij}$ is the antisymmetric tensor, the derivatives
are denoted as $F^{\prime}\equiv\frac{\partial F}{\partial b^{2}_{\perp}}$ and
$\Delta_{b_{\perp}}F\equiv=4\frac{\partial}{\partial b^{2}}\left({\bf
b}_{\perp}^{2}\frac{\partial}{\partial b^{2}}\right)F$ and the GPD
$\bar{E}_{T}=E_{T}+2\widetilde{H}_{T}$. The moments of $\rho$ given by

$\langle x^{n-1}\rangle_{\rho}({\bf b}_{\perp},{\bf s}_{\perp},{\bf
S}_{\perp})\equiv\int_{-1}^{1}\,dx\,x^{n-1}\rho(x,{\bf b}_{\perp},{\bf
s}_{\perp},{\bf S}_{\perp}),$ (2)

are related to generalized form factors (GFFs), e.g.
$A_{Tn0}=\int_{-1}^{1}\,dx\,x^{n-1}H_{T}$,
$B_{Tn0}=\int_{-1}^{1}\,dx\,x^{n-1}E_{T}$,
$\widetilde{A}_{Tn0}=\int_{-1}^{1}\,dx\,x^{n-1}\widetilde{H}_{T}$. The
relevant GFFs are determined from the two lowest Mellin moments by considering
the nucleon matrix elements of the following operators:
$\bullet$ Tensor operator, yielding the first Mellin tensor GFFs

$\displaystyle\langle
N(p^{\prime},s^{\prime})|\bar{q}\sigma^{\mu\nu}q|N(p,s)\rangle=\bar{u}_{N}(p^{\prime},s^{\prime})\left[\sigma^{\mu\nu}\,A^{q}_{T10}(Q^{2})+i\frac{\gamma^{[\mu}\Delta^{\nu]}}{2m_{N}}\right.\,B^{q}_{T10}(Q^{2})$
(3) $\displaystyle+$
$\displaystyle\left.\frac{\overline{P}^{[\mu}\Delta^{\nu]}}{m_{N}^{2}}\,\widetilde{A}^{q}_{T10}(Q^{2})\right]u_{N}(p,s),$

where with $[...]$ we denote antisymmetrization over the indices and
$\bar{P}^{\mu}=\frac{p^{\prime\mu}+p^{\mu}}{2}$. In the forward limit,
$A_{T10}(0)$ gives the tensor charge $g_{T}$, while
$\bar{B}^{q}_{T10}(0)\equiv B^{q}_{T10}(0)+2\widetilde{A}^{q}_{T10}(0)$ is the
anomalous tensor magnetic moment $\kappa_{T}$, which is related to the Boer-
Mulders function $h^{\perp}\sim-\kappa_{T}$.
$\bullet$ Vector one-derivative operator $\mathcal{O}^{\mu\nu}_{\rm
V}=\bar{q}\gamma^{\\{\mu}iD^{\nu\\}}q$, yielding the unpolarized second Mellin
GFFs

$\displaystyle\langle N(p^{\prime},s^{\prime})|\mathcal{O}^{\mu\nu}_{\rm
V}|N(p,s)\rangle=\bar{u}_{N}(p^{\prime},s^{\prime})\Big{[}A_{20}(Q^{2})\gamma^{\\{\mu}\overline{P}^{\nu\\}}+B_{20}(Q^{2})\frac{i\sigma^{\\{\mu\alpha}q_{\alpha}\overline{P}^{\nu\\}}}{2m_{N}}$
(4) $\displaystyle+$ $\displaystyle
C_{20}(Q^{2})\frac{q^{\\{\mu}q^{\nu\\}}}{m_{N}}\Big{]}u_{N}(p,s),$

where with $\\{...\\}$ we denote symmetrization over the indices and
subtraction of the trace. In the forward limit, $A_{20}(0)$ gives the quark
momentum fraction $\langle x\rangle_{q}$ and the total quark spin is
$J_{q}=\frac{1}{2}\left[A_{20}(0)+B_{20}(0)\right]$.
$\bullet$ Tensor one-derivative operator ${\cal
O}_{T}^{\mu\nu\rho}=i\,\bar{q}\sigma^{[\mu\\{\nu]}\,i\,D^{\rho\\}}\,q$,
yielding the second Mellin tensor GFFs with the transversity moment ${\langle
x\rangle_{\delta q}=A^{q}_{T20}(0)}$.

$\displaystyle\langle
N(p^{\prime},s^{\prime})|{\cal{O}}_{T}^{\mu\nu\rho}|N(p,s)\rangle=u_{N}(p^{\prime},s^{\prime})\left[i\sigma^{\mu\nu}\overline{P}^{\rho}\,A_{T20}(Q^{2})+\frac{\gamma^{[\mu}\Delta^{\nu]}}{2m_{N}}\overline{P}^{\rho}\,B_{T20}(Q^{2})\,\right.$
(5) $\displaystyle+$ $\displaystyle
i\,\frac{\overline{P}^{[\mu}\Delta^{\nu]}}{m_{N}^{2}}\overline{P}^{\rho}\,\widetilde{A}_{T20}(Q^{2})+\left.\frac{\gamma^{[\mu}\overline{P}^{\nu]}}{m_{N}}\Delta^{\rho}\,\widetilde{B}_{T21}(Q^{2})\right]u_{N}(p,s).$

##  3 Results

Table 1: The first column gives the ensemble name, the second the lattice volume, the third the lattice spacing, the fourth the value of the pion mass $m_{\pi}$, the fifth the spatial length in fm and the last $m_{\pi}L$. ensemble name |  $(L/a)^{3}.T/a$ |  $a$ (fm) |  $m_{\pi}$ (MeV) |  $L$ (fm) | $m_{\pi}L$
---|---|---|---|---|---
cB211.072.64 | $64^{3}\cdot 128$ | $0.07961~{}(13)$ | $140.2~{}(2)$ | $5.09$ | $3.62$
cC211.060.80 | $80^{3}\cdot 160$ | $0.06821~{}(12)$ | $136.7~{}(2)$ | $5.46$ | $3.78$
cD211.054.96 | $96^{3}\cdot 192$ | $0.05692~{}(10)$ | $140.8~{}(2)$ | $5.46$ | $3.90$

Figure 1: We show the analysis of the ratio from which we extract
$g_{T}^{u-d}$. The left panels show the ratio versus the insertion time
$t_{\rm ins}$ shifted by half the source-sink time $t_{s}/2$. The middle
panels show the values extracted when fitting the ratio on the left with a
constant leaving $8a$ time separation from the source and sink, or the value
of the central point when not enough values of $t_{s}$ are available. In the
right panel, we show the values obtained when we perform a two state-fit (red
circles), a three-state fit (blue triangles) and the summation method (green
circles), as we vary the lower fit range, $t_{s,{\rm min}}$ used in the fits.
The open symbol shows the value we take and the red band in the middle panels
shows the prediction of the plateau values for a given $t_{s}$ value using the
selected parameters determined by the open symbols. Results are shown for the
cB211.072.64 (top), cC211.060.80 (middle) and cD211.054.96 (bottom) ensemble.
For the latter, we are still increasing statistics for the three largest
values of $t_{s}$, see Table. 2.

We illustrate our analysis by considering the extraction of $g_{T}$ for the
three $N_{f}=2+1+1$ twisted mass clover-improved fermion ensembles simulated
by the Extended Twisted Mass Collaboration (ETMC) using physical values of the
two light, the strange and the charm quark masses. The parameters of the
ensembles are listed in Table 1. In Fig. 1, we show the extraction of the
isovector $g_{T}^{u-d}$ using three physical listed in Table 1. As can be
seen, there are large excited state contributions that would require $t_{s}>2$
fm to eliminate in the ratio [4]. However, including a the first excited
state, the results converge rapidly yielding the nucleon matrix element. This
is confirmed by including the second excited state in the fit. Note that
having approximately the same statistical error as we increase $t_{s}$ is
important for identifying clearly the nucleon matrix element. The summation
method, on the other hand, shows slow convergence, but eventually converges to
the chosen value from the two-state fit. In Table 2 we give the statistics
used.

cB211.072.64
---
750 configurations
$t_{s}/a$ |  $t_{s}$[fm] | $n_{src}$
8 | 0.64 | 1
10 | 0.80 | 2
12 | 0.96 | 4
14 | 1.12 | 6
16 | 1.28 | 16
18 | 1.44 | 48
20 | 1.60 | 64
Nucleon 2pt | 264

cC211.060.80
---
400 configurations
$t_{s}/a$ |  $t_{s}$[fm] | $n_{src}$
6 | 0.41 | 1
8 | 0.55 | 2
10 | 0.69 | 4
12 | 0.82 | 10
14 | 0.96 | 22
16 | 1.10 | 48
18 | 1.24 | 45
20 | 1.37 | 116
22 | 1.51 | 246
Nucleon 2pt | 650

cD211.054.96
---
500 configurations
$t_{s}/a$ |  $t_{s}$[fm] | $n_{src}$
8 | 0.46 | 1
10 | 0.57 | 2
12 | 0.68 | 4
14 | 0.80 | 8
16 | 0.91 | 16
18 | 1.03 | 32
20 | 1.14 | 64
22 | 1.25 | 16
24 | 1.37 | 32
26 | 1.48 | 64
Nucleon 2pt | 368

Table 2: Statistics for computing the isovector matrix element of the tensor
charge for the cB211.072.64 (top), the cC211.060.80 (middle) and the
cD211.054.96 (bottom) ensemble. The third column gives the number of source
positions used. For the cD211.054.96 ensemble, we are increasing statistics
for $t_{s}/a=22,\,24$ and 26.

Using the results for the three physical mass point ensembles we extrapolate
to the continuum limit avoiding any chiral extrapolations that may introduce
unquantified systematic errors. We obtain a value of $g_{T}^{u-d}=0.924(54)$
[5]. In Fig. 2, we compare our result with those from other lattice QCD
collaborations, all of which have used ensembles simulated with larger than
physical pion mass and chirally extrapolated the results to the physical
point. As can be seen, the lattice QCD results are in agreement and more
accurate than the values extracted from phenomenology, demonstrating that
lattice QCD can provide a precise determination.

Figure 2: We show lattice QCD results (green symbols for $N_{f}=2+1+1$ and
yellow for $N_{f}=2+1$ ensembles) for $g_{T}^{u-d}$ from ETMC [5, 6] using
only physical mass point ensembles; PNDME [7]; $\chi$QCD [8]; NME [9]; and
CLS-Mainz [10]. Open green and yellow symbols denote lattice QCD results
without continuum extrapolation. Results from phenomenology are shown with the
black circles [11, 12, 13, 14].

For the cB211.072.64 ensemble, we also calculate the tensor charge for the
each quark flavor. This requires the computation of disconnected
contributions. We use hierachical probing to improve the signal and 600
thousand two-point functions. For the light quark loops we use, in addition,
deflation of the 200 lowest modes. The resulting values are given in Table 3.

| ${u\texttt{-}d}$ | ${u}$ | ${d}$ | $s$ | $c$
---|---|---|---|---|---
$g_{T}$ | 0.936(25) | 0.729(22) | -0.207(75) | -0.00268(58) | -0.00024(16)
Table 3: We give our values for $g_{T}$ extracted using the cB211.072.64
ensemble: first column for $g^{u-d}_{T}$ and the 2nd to 5th for $g^{u}_{T}$,
$g_{T}^{d}$, $g_{T}^{s}$, $g_{T}^{c}$, respectively [6]. The values are given
in the $\overline{\rm MS}$ scheme at 2 GeV.

The second Mellin GFFs are extracted by performing a similar analysis. As an
example, we show in Fig. 3 the extraction of the isovector tensor moment
$\langle x\rangle_{\delta u-\delta d}=A^{u-d}_{T20}(0)$ for the cC211.060.80
ensemble.

Figure 3: Excited states analysis of the ratio to extract $A_{T20}$ for the
cC211.060.80 ensemble. The notation is the same as that of Fig. 1

Having results for the three ensembles with different lattice spacings, we can
perform a continuum extrapolation linearly in $a^{2}$. For the GFFs the
continuum extrapolation is done at a fixed value of $Q^{2}$ interpolating the
$Q^{2}$ of two of the ensembles to match the third one. The continuum
extrapolations are shown in Fig. 4.

Figure 4: Left: We show the continuum extrapolation of the isovector GFFs for
$Q^{2}=0$. Right: We show the fits of the $Q^{2}$-dependence of the GFFs after
extrapolated to the continuum limit. $A_{T10}(0)=g^{u-d}_{T}$, $A_{20}(0)$ is
the isovector momentum fraction and $A_{T20}=\langle x\rangle_{\delta u-\delta
d}$. Figure 5: Contours of $\int_{-1}^{1}\rho(x,{\bf b}_{\perp},{\bf
s}_{\perp},{\bf S}_{\perp})=\langle 1\rangle_{\rho}$ [fm-2] as a function of
$b_{x}$ and $b_{y}$. Top-left: transversely polarized quarks in an unpolarized
nucleon; top-right: unpolarized quarks in a transversely polarized nucleon;
bottom-left: transversely polarized quarks in a transversely polarized
nucleon; and bottom-right: transversely polarized quarks in a perpendicularly
polarized nucleon.  Figure 6: We show $\int_{-1}^{1}\,dxx\rho(x,{\bf
b}_{\perp},{\bf s}_{\perp},{\bf S}_{\perp})=\langle x\rangle_{\rho}$ [fm-2].
The notation is the same as in Fig. 5.

We fit the $Q^{2}$-dependence of the isovector GFFs after taking the continuum
limit using the form $F(Q^{2})=\frac{F(O)}{1+Q^{2}/m^{2})^{p}}$, where $F(0)$,
$m$ and $p$ are fit parameters. We Fourier transform into impact parameter
with the function [3]

$F(b^{2}_{\perp})=\frac{m^{2}F(0)}{2^{p}\pi\Gamma(p)}(mb_{\perp})^{p-1}K_{p-1}(mb_{\perp}),$
(6)

where $\Gamma(x)$ is the Euler gamma function, $K_{n}(x)=K_{-n}(x)$ the
modified Bessel functions and $b_{\perp}=\sqrt{b_{\perp}^{2}}$. This then
allows us to compute the first and second moments of the transverse density
distributions, shown in Figs. 5 and 6.

We observe a sizeable deformation for the first moment. We consider four
cases: i) Transversely polarized quarks in an unpolarized nucleon, where we
observe a mild distortion because in the isovector combination the
$\bar{B}_{T10}(b^{2}_{\perp})$ term has a milder $Q^{2}$-dependence; ii)
Unpolarized quarks in a transversely polarized nucleon, where we observe a
large distortion towards the positive $b_{y}$ direction. This can be traced
back to the GFF $B_{10}$, contributing to the term for $E^{\prime}$ that, as
seen in Fig. 4, is rapidly decreasing resulting in a large distortion. The
origin of this behavior is related to the Sivers effect [16], a connection
that was made in Ref. [15]; iii) Transversely polarized quarks in a
transversely polarized nucleon, where all the terms in Eq. (1) contribute
leading to a significant deformation of the density; and iv) Transversely
polarized quarks in a perpendicularly polarized nucleon, where the third term
in Eq. (1) drops out and the fourth one creates a significant impact, leading
to a large distortion in the $b_{x}$ direction. For the second moment $\langle
x\rangle_{\rho}$ the distortion in all cases is reduced.

##  4 Conclusions

We show that the moments of PDFs can be extracted precisely and directly at
the physical pion mass. Taking the continuum limit at fixed physical pion mass
removes any uncontrolled errors due to the chiral extrapolation. The value of
the isovector tensor charge of the nucleon provides valuable input to
phenomenological analysis, as demonstrated e.g. by the the recent analysis by
the JAM collaboration [17], which with lattice input can improve their
determination of transverse spin asymmetries. The calculation of sea quark
contributions is feasible using physical pion mass ensembles providing, as
demonstrated in this work, the tensor charge for each quark flavor. Although
the disconnected contributions are small, they are non-zero for the light,
strange and charm quark flavors. Using the computation of the unpolarized and
transversity GFFs in the continuum limit, one can determine the moments of the
transverse quark spin density. We find large deformation for the first moment
showing a sizeable Sivers effect. We determine the anomalous magnetic moment
$\kappa_{T}\equiv\bar{B}_{T10}(0)=1.051(94)$, a fundamental quantity that
describe the deformation of the transverse polarized quark distribution in an
unpolarized nucleon. First lattice results on $\kappa_{T}$ were presented by
the QCDSF/UKQCD collaboration [18]. Since $\kappa_{T}\sim-h_{1}^{\perp}$ [19],
our value confirms that the Boer-Mulders function, $h_{1}^{\perp}$, is
negative and sizeable and corroborates the result found in the lattice QCD
study of Ref. [20] for the transverse momentum dependent PDFs.

## Acknowledgments

We would like to thank all members of the Extended Twisted Mass Collaboration
for a very constructive and enjoyable collaboration. M.C. acknowledges
financial support by the U.S. Department of Energy Early Career Award under
Grant No. DE-SC0020405. K.H. is supported by the Cyprus Research and
Innovation foundation under contract number POST-DOC/0718/0100 and CULTURE-
AWARD-YR/0220/0012. S.B., J.F., K.H. and G. S. are partly supported by the
EuroCC project (GA No. 951732) and by the H2020 project PRACE 6-IP (GA No.
823767). The project acknowledges support from the European Joint Doctorate
projects HPC-LEAP and STIMULATE funded by the European Union’s Horizon 2020
research and innovation programme under grant agreement No 642069 and 765048,
respectively. P.D. acknowledges support from the European Union’s Horizon 2020
research and innovation programme under the Marie Sklodowska-Curie grant
agreement No. 813942 (EuroPLEx) and from INFN under the research project INFN-
QCDLAT. The authors gratefully acknowledge computing time allocated on Piz
Daint by CSCS via the projects with ids s702, s954 and pr79 and on JUWELS at
Jülich Supercomputing Centre (JSC) through the projects FSSH, PR74YO and
CECY00, CHCH02. Part of the results have been produced within the EA program
of JUWELS Booster also with the help of the JUWELS Booster Project Team (JSC,
Atos, ParTec, NVIDIA). We acknowledge the CINECA award under the ISCRA
initiative, for the availability of high performance computing resources and
support and PRACE for awarding us access to Marconi100 at CINECA (Italy), Piz-
Daint at CSCS (Switzerland) and Hawk at HLRS (Germany). We acknowledge the
Gauss Centre for Supercomputing e.V. for project pr74yo providing computing
time on SuperMUC at LRZ and Juwels Booster at JSC. The authors acknowledge the
Texas Advanced Computing Center (TACC) at The University of Texas at Austin
for providing HPC resources (Project ID PHY21001).

## References

  * [1] A. Accardi, J. L. Albacete, M. Anselmino, N. Armesto, E. C. Aschenauer, A. Bacchetta, D. Boer, W. K. Brooks, T. Burton and N. B. Chang, et al. Eur. Phys. J. A 52, no.9, 268 (2016) arXiv:1212.1701.
  * [2] G. Martinelli and C. T. Sachrajda, Phys. Lett. B 217, 319-324 (1989).
  * [3] M. Diehl and P. Hagler, Eur. Phys. J. C 44, 87-101 (2005) arXiv:hep-ph/0504175.
  * [4] O. Bär, Phys. Rev. D 94, no.5, 054505 (2016) arXiv:1606.09385.
  * [5] C. Alexandrou, S. Bacchio, M. Constantinou, P. Dimopoulos, J. Finkenrath, R. Frezzotti, K. Hadjiyiannakou, K. Jansen, B. Kostrzewa and G. Koutsou, et al. arXiv:2202.09871.
  * [6] C. Alexandrou, S. Bacchio, M. Constantinou, J. Finkenrath, K. Hadjiyiannakou, K. Jansen, G. Koutsou and A. Vaquero Aviles-Casco, Phys. Rev. D 102, no.5, 054517 (2020) arXiv:1909.00485.
  * [7] R. Gupta, Y. C. Jang, B. Yoon, H. W. Lin, V. Cirigliano and T. Bhattacharya, Phys. Rev. D 98, 034503 (2018) arXiv:1806.09006.
  * [8] D. Horkel et al. [$\chi$QCD], Phys. Rev. D 101, no.9, 094501 (2020) arXiv:2002.06699.
  * [9] S. Park et al. [Nucleon Matrix Elements (NME)], Phys. Rev. D 105, no.5, 054505 (2022) arXiv:2103.05599.
  * [10] T. Harris, G. von Hippel, P. Junnarkar, H. B. Meyer, K. Ottnad, J. Wilhelm, H. Wittig and L. Wrang, Phys. Rev. D 100, no.3, 034513 (2019), arXiv:1905.01291.
  * [11] M. Radici and A. Bacchetta, Phys. Rev. Lett. 120, no.19, 192001 (2018) arXiv:1802.05212.
  * [12] Z. B. Kang, A. Prokudin, P. Sun and F. Yuan, Phys. Rev. D 93, no.1, 014009 (2016) arXiv:1505.05589.
  * [13] G. R. Goldstein, J. O. Gonzalez Hernandez and S. Liuti, arXiv:1401.0438. [14]
  * [14] M. Pitschmann, C. Y. Seng, C. D. Roberts and S. M. Schmidt, Phys. Rev. D 91, 074004 (2015) arXiv:1411.2052.
  * [15] M. Burkardt, Nucl. Phys. A 735, 185-199 (2004) arXiv:hep-ph/0302144.
  * [16] M. Bury, A. Prokudin and A. Vladimirov, JHEP 05, 151 (2021) arXiv:2103.03270.
  * [17] L. Gamberg et al. [Jefferson Lab Angular Momentum (JAM) and Jefferson Lab Angular Momentum], Phys. Rev. D 106, no.3, 034014 (2022) arXiv:2205.00999.
  * [18] M. Göckeler et al. [QCDSF and UKQCD], Phys. Rev. Lett. 98, 222001 (2007) arXiv:hep-lat/0612032.
  * [19] M. Burkardt, Phys. Rev. D 72, 094020 (2005) arXiv:hep-ph/0505189.
  * [20] B. U. Musch, P. Hagler, M. Engelhardt, J. W. Negele and A. Schafer, Phys. Rev. D 85, 094510 (2012) arXiv:1111.4249.


# The Persian Dependency Treebank Made Universal

Mohammad Sadegh Rasooli1 Rasooli and Safari equally contributed in the
conversion and experimentation process. Rasooli and Moloodi equally
contributed in the linguistic design of conversion rules and manual
investigation of conversions.  Pegah Safari2* Amirsaeid Moloodi3* Alireza
Nourian4
1 Department of Computer and Information Science, University of Pennsylvania,
Philadelphia, PA, USA
2 Faculty of Computer Science and Engineering, Shahid Beheshti University,
Tehran, Iran
3 Department of Foreign Languages and Linguistics, Shiraz University, Shiraz,
Iran
4 Sobhe, Tehran, Iran
<EMAIL_ADDRESS><EMAIL_ADDRESS>
<EMAIL_ADDRESS><EMAIL_ADDRESS>

###### Abstract

We describe an automatic method for converting the Persian Dependency Treebank
Rasooli et al. (2013) to Universal Dependencies. This treebank contains
$29107$ sentences. Our experiments along with manual linguistic analysis show
that our data is more compatible with Universal Dependencies than the Uppsala
Persian Universal Dependency Treebank Seraji et al. (2016), and is larger in
size and more diverse in vocabulary. Our treebank brings in a labeled
attachment F-score of $85.2$ in supervised parsing. Our delexicalized Persian-
to-English parser transfer experiments show that a parsing model trained on
our data is $\approx$2% absolutely more accurate than that of Seraji et al.
(2016) in terms of labeled attachment score.111The treebank is publicly
available in https://github.com/UniversalDependencies/UD_Persian-
PerDT/tree/dev.

##  1 Introduction

In recent years, there has been a great deal of interest in developing
universal dependency treebanks McDonald et al. (2013); Rosa et al. (2014);
Nivre et al. (2020). The main goal of the Universal Dependencies project Nivre
et al. (2020) is to develop a consistent linguistic annotation scheme in
different levels from tokenization to syntactic dependency relations. As a
result, the majority of annotation discrepancies disappear, and the resulting
dataset facilitates several cross-lingual natural language processing tasks
including part-of-speech transfer Täckström et al. (2013), syntactic transfer
Naseem et al. (2010); McDonald et al. (2011); Ammar et al. (2016); Zhang et
al. (2019), and probing Tenney et al. (2019); Hewitt and Manning (2019).
Starting with 10 treebanks in 2015, there are 163 treebanks in version 2.6
(May 2020) including the Uppsala Persian Treebank Seraji et al. (2016).

Persian (aka Farsi) is a pro-drop morphologically rich language with a high
degree of free word order and a unique light verb construction Karimi-Doostan
(2011). Despite its importance, it still suffers from lack of sufficient
annotated data. The Uppsala Universal treebank Seraji et al. (2016) is
currently the only publicly available universal treebank for Persian. It is a
valuable resource based on news genre, and has been used as a testbed in
previous work Zeman et al. (2018); Chi et al. (2020). Among other non-
universal treebanks, the Persian dependency treebank (PerDT) Rasooli et al.
(2013) is significantly larger than Seraji et al. (2016) (29K vs. 6K
sentences), and its sentences are sampled from contemporary Persian texts in
different genres (as opposed to only news genre).

Figure 1: An example of our automatic conversion. The universal labels are
shown at the bottom of words.

In this paper, we propose an automatic method for converting PerDT Rasooli et
al. (2013) to Universal Dependencies (An example of such conversion is shown
in Figure 1). After a thorough analysis of dependency relations in the
treebank, we design different mapping rules to generate trees with universal
relations. This process involves a series of steps including unifying
tokenization, part-of-speech tags, named-entity recognition, and finally
mapping dependencies. The mapping for many relations are not necessarily one-
to-one, and we have to deal with peculiar cases that are specific to certain
structures in modern Persian. Therefore, our approach is neither a blind one-
to-one mapping, nor an expensive and time-consuming manual process. We
empirically show that our annotations are more compatible with the Universal
guidelines via learning a delexicalized transfer model with more than $2\%$
absolute difference in labeled attachment score. The summary of our
contributions is as following:

  * •

We propose an automatic annotation conversion process with manual care of
special cases. We develop a new Persian Universal Treebank with 29107
sentences. This is in contrast to the treebank of Seraji et al. (2016) that
contains 5997 sentences.

  * •

We develop a modified and corrected version of PerDT with the Universal
tokenization scheme. Moreover, the new release resolves various tagging errors
in the original dataset. Most of these corrections are made by manually fixing
annotation errors flagged by our mapping pipeline.

##  2 Related Work

There has been a great deal of interest in designing and developing Persian
dependency treebanks Pouramini and Mozayani (2007); Seraji et al. (2012, 2014,
2016); Rasooli et al. (2011b, 2013); Ghayoomi and Kuhn (2014). Among them, the
Uppsala UD treebank Seraji et al. (2016) is the only treebank with Universal
Dependencies. We have found some caveats in the Uppsala Universal Treebank
Seraji et al. (2016). This causes annotation discrepancies in some frequently
used dependency relations such as _compound:lvc_ , _cop_ , _csubj_ , _fixed_ ,
_obl_ , and _xcomp_ (see §A for more details).

We primarily focus on converting the Persian dependency treebank (PerDT)
Rasooli et al. (2013). PerDT has been used in previous studies for Persian
dependency parsing Khallash et al. (2013); Feely et al. (2014); Nourian et al.
(2015); Pakzad and Minaei-Bidgoli (2016). It has been extended to other
representations including semantic roles Mirzaei and Moloodi (2016) and
discourse Mirzaei and Safari (2018). It is also included in the HamleDT
collection Rosa et al. (2014).

##  3 Approach

In the conversion process, we have noticed several key differences between
PerDT and UD. We decompose the conversion process into 3 steps: 1)
tokenization, 2) part-of-speech mapping, 3) systematic changes to PerDT, and
4) dependency relation mapping. In this section, we briefly describe the
mentioned steps.

###  3.1 Tokenization

There are two key differences in PerDT tokenization from UD: 1) Multiword
inflections of simple verbs in Persian are grouped as one word with spaces in
between parts following the deterministic rules from Rasooli et al. (2011a).
We follow the guidelines in (Rasooli et al., 2013, Table 3) to find out the
main verb and make other parts an “aux” dependent of the main verb. We
introduce the “AUX” part-of-speech tag and “aux” dependency relation
(“aux:pass” for passive verbs) in this UD tokenization scheme. 2) Clitics are
only detached from words in cases for which they play an object or verbal
role. Other clitics are pronominal clitics attached to nouns, prepositions,
pronouns and adjectives. By looking at the word lemma, we recover those
pronouns, and detach them, and assign their heads to the closest nominal word
with the “MOZ” (_Ezafe_) dependency label.

###  3.2 POS Mapping

This is the most straightforward step except for proper nouns. We could only
discover a small portion of them by finding noun phrases with an identifier
(IDEN POS tag for words such as “Dr.” or “Mr."). In addition to mapping the
IDEN POS to PROPN, we use a recent BERT-based Persian named-entity tagger
Taher et al. (2020) to recover additional proper nouns. The tagger can find 7
different entities including date, location, money, organization, percent,
person and time. We only consider the _person_ and _location_ entities, and
manually revise the results to add missing entities, foreign words and the
name of months. Table 1 shows the mappings.

PerDT | Condition | UD
---|---|---
V |  | VERB
N | NER=False | NOUN
NER=True | PROPN
SUBR |  | SCONJ
CONJ |  | CCONJ
ADV |  | ADV
ADJ | NER=False | ADJ
NER=True | PROPN
PR |  | PRON
PUNC |  | PUNCT
ADR |  | INTJ
IDEN |  | PROPN
PART | Word= rA > | ADP
Word$\in${xwb>, Axr>} | INTJ
Otherwise | PART
PREM |  | DET
PRENUM | Cardinal | NUM
Ordinal | ADJ
PREP |  | ADP
POSTP |
POSTNUM | Cardinal | NUM
Ordinal | ADJ
PSUS |  | INTJ
Table 1: Mapping rules for part-of-speech tags.  PerDT | Precondition | Pre-action | UD
---|---|---|---
ACL |  | CMR (mark) [fig.3(b)] | ccomp
PRD |
ADV | ?$\rightarrow${ADJ,ADV} |  | advmod
Otherwise |  | obl
AJCONJ |  | Conj rotation [fig.3(a)] | conj
AVCONJ |
NCONJ |
PCONJ |
VCONJ |
AJUCL |  | CMR (mark) [fig.3(b)] | advcl
| ?$\rightarrow${ADP} | CMR (case) [fig.3(b)] | obl
| ?$\rightarrow${NUM} |  | nummod
APOSTMOD | ?$\rightarrow$_nominal_ |  | nmod
APREMOD | ?$\rightarrow${ADJ} |  | amod
| ?$\rightarrow$DET |  | det
| ?$\rightarrow${ADV} |  | advmod
ADVC |  | CMR (case) [fig.3(b)] | obl:arg
NEZ |
AJPP |
VPP |
VPRT |
APP |  |  | appos
COMPPP |  $\exists$ dep | CMR (case) [fig.3(b)] | case
Otherwise |  | fixed
NEZ |  |  |
ENC |  |  | compound:lvc
NE |  |
NPRT |  | CMR (case) [fig.3(b)]
NVE |  |
LVP |  |  | compound:lv
NCL | ?$\rightarrow$SCONJ | CMR (mark) [fig.3(b)] | acl
Otherwise |
MESU |  | Dep $\rightarrow$Head (Flip) | nmod
MOS | AUX$\rightarrow$? | Dep $\rightarrow$Head (Flip) | cop
Otherwise |  | xcomp
MOZ | ?$\rightarrow$NOUN | CMR (case) [fig.3(b)] | nmod
?$\rightarrow$ADJ | amod
Otherwise | advmod
NADV | ?$\rightarrow$NOUN | CMR (case) [fig.3(b)] | nmod
?$\rightarrow$ADJ | amod
Otherwise | advmod
NPOSTMOD |  |  | amod
NPP | {NVE$|$ENC}$\rightarrow$? | NPP rotation [fig.3(c)] | obl:arg
Otherwise | CMR (case) [fig.3(b)] | nmod
NPREMOD | ?$\rightarrow$DET |  | det
?$\rightarrow$ cardinal |  | nummod
Otherwise |  | amod
OBJ |  $\exists$ OBJ2 sib. |  | iobj
Otherwise |  | obj
OBJ2 |  |  | obj
PARCL | ?$\rightarrow$CCONJ | Dep $\rightarrow$Head (Flip) | conj
Otherwise | parataxis
PART |  |  | mark
PUNC |  |  | punct
PROG | Active$\rightarrow$? |  | aux
Passive$\rightarrow$? |  | aux:pass
ROOT |  |  | root
SBJ | Active$\rightarrow$? |  | nsubj
Passive$\rightarrow$? |  | nsubj:pass
TAM |  |  | xcomp
VCL | Modal verb$\rightarrow$? | Dep $\rightarrow$Head (Flip) | aux
$\exists$ MOS,$\nexists$ SUBJ sib. | CMR (mark) [fig.3(b)] | csubj
Otherwise | ccomp
PREDEP | NUM$\rightarrow$? |  | advmod
NOUN$\rightarrow$PRON |  | dislocated
?$\rightarrow$CCONJ |  | cc
?$\rightarrow$NOUN |  | obl
_Last mapping_ |  | advmod
POSDEP | NOUN$\rightarrow${hm, nyz>} |  | dep
?$\rightarrow$NOUN |  | obl
?$\rightarrow$CCONJ |  | cc
_Last mapping_ |  | advmod
Table 2: Mapping rules for dependencies. PerDT labels are described in Rasooli
et al. (2013, Table 2). First _Preconditions_ (2nd column) should satisfy.
Afterwards, _Preactions_ (3rd column) are applied before applying the UD
conversions (4th column). These preactions are depicted in Figure 3.

###  3.3 Systematic Changes to PerDT

Before starting to convert the treebank Dadegan Research Group (2012), we have
made the following systematic changes to PerDT:

Figure 2: The result of applying rotations of conjunctions for the sentence
_Why didn’t you defend in the court and acted like that?_. In this example,
case rotation for preposition is also shown.

  * •

We convert the order of verbal conjunctions in the original data. In PerDT,
verbal conjunctions are conventionally attached from the end to the beginning
Dadegan Research Group (2012).222Examples in https://bit.ly/2Mfz1iH We find
this convention unintuitive and reverted the order of conjunctions. Figure 2
shows an example of such rotation.

  * •

Words such as “billion”, “million”, “thousand” are tagged as nouns.333Examples
in https://bit.ly/2Y105Yv This might be due to the fact that these words can
be inflected as plurals while number should not be inflected in Persian. We
believe that a better tagging decision for these words is number since their
inflection as plurals is due to a special kind of _zero derivation_ or
_conversion_ numbers to nouns in particular contexts Booij (2012).

  * •

PerDT assumes that all inflections of “^sdn>” [Sodæn] is passive and its lemma
is “krdn>” [kædæn]. We have changed this assumption and use the superficial
lemma for those instances. The decision makes our data similar to the
annotations of Seraji et al. (2016).

Correction Type | # | %
---|---|---
| Systematic | 3694 | 0.762
Lemma | Others | 59 | 0.012
| Systematic | 529 | 0.109
POS | Others | 298 | 0.061
| Systematic | 3693 | 0.762
FPOS | Others | 90 | 0.018
| Systematic | 27407 | 5.658
Dependency head | Others | 967 | 0.199
| Systematic | 18516 | 3.823
Dependency label | Others | 656 | 0.135
Word Form | 39 | 0.008
Table 3: Statistics of PerDT corrections. By systematic, we mean deterministic
corrections such as verbal conjunctions (see §3.3 for details).

Table 3 shows the statistics of changes that we have made to the data
including systematic changes and fixes to incorrect annotations.

###  3.4 Dependency Relation Mapping

PerDT contains 43 syntactic relations for which many of them cannot easily map
to UD. Moreover, conjunctions in PerDT are arranged from the beginning of the
sentence to the end in chain-style manner. More importantly, compared to UD
scheme for which content words are considered as heads, PerDT assigns
prepositions as the head of prepositional phrases and auxiliary verbs as the
head of sentences.

Before applying the conversion rules, we label words that are not well-edited
and typed as more than one token as _goeswith_. We then label proper noun
phrases that are not syntactically compositional as _flat:name_. We also
analyze complex numbers as _flat:num_ and their coordinating conjunctions as
_cc_ dependent of each following word. Afterwards we follow the rules in Table
2. As depicted in the Table, there are conditions that should be satisfied
before applying a conversion, and some actions such as flipping a head with
its dependent are needed before certain mappings. Finally, we label the few
remaining undecided dependencies as _dep_.

(a) Conj Rotation

(b) CMR (role)

(c) NPP rotation

Figure 3: A graphical depiction of rotation rules used in this work (see Table
2 for their use cases).

##  4 Experiments and Analysis

The general statistics of our data vs. the Uppsala treebank Seraji et al.
(2016) are shown in Table 4. We observe that our data is superior in many
aspects including size and diversity compared to the Uppsala Treebank Seraji
et al. (2016). The most important fact about PerDT is that its sentences are
intentionally sampled in order to cover almost all verbs from the Verb Valency
Lexicon Rasooli et al. (2011b) leading to 3.9 times more verb lemmas than the
Uppsala Treebank. Table 5 shows the counts of each dependency label in the
converted Data.

Part | Data | Sen.# | Tok.# | Types#
---|---|---|---|---
Word | Lemma | Verb
Train | UDT | 4798 | 122K | 13.9K | 6.7K | 1226
Ours | 26196 | 459K | 34.9K | 20.7K | 5275
Dev | UDT | 599 | 15K | 3.9K | 2.0K | 278
Ours | 1456 | 26K | 7.0K | 5.2K | 1427
Test | UDT | 600 | 16K | 3.9K | 3.1K | 385
Ours | 1455 | 24K | 6.7K | 5.1K | 1671
All | UDT | 5997 | 154K | 15.8K | 7.6K | 1387
Ours | 29107 | 509K | 36.7K | 21.6K | 5413
Table 4: Statistics of our data vs. UDT Seraji et al. (2016) in different data splits.  Label | Frequency | %
---|---|---
case | 71118 | 14.1
conj | 23739 | 4.7
acl | 10034 | 1.9
obl | 30737 | 6.1
punct | 44336 | 8.8
cop | 6366 | 1.2
det | 10273 | 2
advmod | 9158 | 1.8
aux:pass | 822 | 0.1
nmod | 59442 | 11.6
appos | 1059 | 0.2
aux | 12886 | 0.16
amod | 22576 | 4.4
compound:lvc | 32339 | 6.4
nsubj:pass | 822 | 0.1
nsubj | 27181 | 5.4
name:flat | 7899 | 1.5
dep | 2035 | 0.4
cc | 21300 | 4.2
root | 29107 | 5.8
advcl | 4228 | 0.8
obj | 19999 | 3.9
xcomp | 4920 | 0.9
parataxis | 82 | 0.01
ccomp | 6945 | 1.3
obl:arg | 21510 | 4.2
flat:num | 607 | 0.1
nummod | 5459 | 1
mark | 11982 | 2.3
fixed | 144 | 0.02
compound:lv | 439 | 0.08
csubj | 682 | 0.1
vocative | 174 | 0.03
compound | 42 | 0.008
iobj | 6 | 0.001
dislocated | 1 | 0.0001
Table 5: Frequency of each universal label the converted dataset.

#### Supervised Parsing

We evaluate the resulting data by training UDPipe V.2 Straka and Straková
(2017) along with the pre-trained fastText Grave et al. (2018) embeddings on
our data. We also evaluate our models on the Uppsala treebank Seraji et al.
(2016). Table 6 shows the parsing results using a trained model on our data
and the Uppsala Treebank evaluated by the CoNLL 2018 shared task evaluation
scripts Zeman et al. (2018). It is worth noting that the goal of this
evaluation is not to show which dataset brings in better parsing accuracy: it
is clear that the bigger the dataset is, the higher the accuracy can be. Our
goal is to show that there is a significant performance difference between the
models trained on the two datasets by using the exact same training pipeline.
As shown in Table 6, we see that there is a huge tagging and parsing
performance difference when we move across the datasets. There are two
possible reasons: domain mismatch, and annotation discrepancy. Our analysis
show that annotation discrepancy plays an important role here. As described in
§A, there are some core incompatibilities between the Uppsala treebank Seraji
et al. (2016) and Universal Dependencies guidelines. Our detailed analysis
shows that most of cross-dataset errors come from errors in _nmod_ , _obl_ ,
_fixed_ , and _xcomp_. This is in fact consistent with our manual analysis in
§A.

Test Data | PerDT (Ours) | Seraji et al. (2016)
---|---|---
ID tagger | $\times$ | ✓ | $\times$ | ✓
ID parser | $\times$ | $\times$ | ✓ | $\times$ | $\times$ | ✓
Tokens | 99.9 | 99.99 | 100 | 100.0
Words | 99.1 | 99.64 | 99.7 | 99.59
UPOS | 82.9 | 96.11 | 81.9 | 95.75
Lemmas | 80.7 | 96.20 | 90.2 | 89.55
UAS | 71.2 | 71.2 | 88.4 | 69.5 | 69.8 | 83.5
LAS | 64.4 | 62.6 | 85.2 | 62.1 | 61.0 | 79.4
CLAS | 59.9 | 59.3 | 81.6 | 56.9 | 56.1 | 74.8
MLAS | 49.5 | 54.5 | 78.9 | 46.0 | 53.9 | 73.0
BLEX | 44.6 | 56.9 | 78.2 | 52.1 | 49.0 | 65.5
Table 6: Parsing results based on the CoNLL shared task 2018 Zeman et al.
(2018) evaluation. ID stands for in-domain for which the same training set is
used for training a UDPipe model Straka and Straková (2017)

#### Delexicalized Model Transfer

One way to verify our claim about increased consistency of our UD conversion
with the UD guidelines is to learn a transfer model. In this setting, we
follow the delexicalized parser transfer approach which have been extensively
used in previous work Zeman and Resnik (2008); McDonald et al. (2011);
Täckström et al. (2012). We sample the same number of tokens as of Seraji et
al. (2016) from PerDT. Afterwards, we delexicalize both of the treebanks, and
learn a parser using the Yara Parser Rasooli and Tetreault (2015). We train
two models with 15 epochs and evaluate them on the delexicaled test set of the
Universal English Web Treebank Silveira et al. (2014). The model trained on
PerDT significantly outperforms the other model by 2% both in unlabeled and
labeled attachment score (47.31 vs 45.37 UAS, 38.59 vs. 36.45 LAS). This is a
strong indicator that our data is more compatible with the UD annotations.

##  5 Conclusion

We have introduced our approach in making PerDT Rasooli et al. (2013)
universal. During this process, we have faced different challenges such as
annotation errors in the original data, tokenization inconsistencies, lack of
named entities, part-of-speech and dependency label mapping. Due to automatic
conversions and potential annotation errors in the original treebank, there is
always a chance of some annotation incompatibilities between our treebank and
the Universal guidelines. Therefore, we cannot claim that our conversion is
perfect. However, our experiments have shown that our data is more compatible
with the Universal Dependencies guidelines than the Uppsala treebank Seraji et
al. (2016).

## Acknowledgments

We would like to thank Morteza Rezaei-Sharifabadi for his efforts in getting
copyright permissions and making our treebank publicly available. We also
thank Daniel Zeman for helping us in the process of releasing this treebank.

## Appendix

##  Appendix A Problems in the Universal Annotations of the Uppsala Universal
Treebank

We briefly mention some of the problems in the Uppsala Universal Treebank
Seraji et al. (2016):

  * •

Seraji et al. (2016) does not determine the _csubj_ label in their analysis.
For example, in “lAzem Pæst Pu beresæd” (it is necessary for him to arrive),
it is obvious that what comes after “Pæst” is the clausal subject of the
adjectival sentence predicate “lAzem”. A simple syntactic test supports this
viewpoint: one can convert the clausal complement “Pu beresæd” to a noun
phrase “residæn-e u” (his arrival). The new phrase plays the _nsubj_ role of
the sentence. Therefore, the clausal complement of the sentence should be
_csubj_. Our converted data contains $682$ cases of _csubj_.

  * •

Seraji et al. (2016) considers prepositional and possessive complements of
adjectival heads as _nmod_ and _nmod:poss_ respectively. Their analysis
clearly stands in contradiction to UD annotation guideline in which _nmod_ is
used just for dependents of a nominal head. _obl_ is much better suited for
these cases.

  * •

Seraji et al. (2016) consider “‘^sdn>” [Sodæn] (to become) as copula. What UD
asserts under the _cop_ (copula) label is that _“the equivalents of to become
are not copulas despite the fact that traditional grammar may label them as
such.”_ Instead, it should be deemed as a verbal predicate and its second
complement as _xcomp_.

  * •

“ pydA krdn>” (“peydA kærdæn”) and “ .hA.sl krdn>” (“hAsel kærdæn”) are
considered as two-word light verbs Moloodi and Kouhestani (2017). We consider
the non-verbal part as the first part of the two-word light verb, and use the
_compound:lv_ label for it ($439$ cases in PerDT). However, Seraji et al.
(2016) annotate the nonverbal elements of these complex predicates as _obj_
and considers “peydA” as a nonverbal element.

  * •

_iobj_ label is absent in Seraji et al. (2016), most likely due to the low
frequency of this syntactic relation. Our converted treebank contains $6$
cases of _iobj_.

  * •

Proper nouns are not labeled in Seraji et al. (2016). Ours covers proper nouns
(more than $23K$ tokens).

## References

  * Ammar et al. (2016) Waleed Ammar, George Mulcaire, Miguel Ballesteros, Chris Dyer, and Noah A. Smith. 2016.  Many languages, one parser.  _Transactions of the Association for Computational Linguistics_ , 4:431–444.
  * Booij (2012) Geert Booij. 2012.  _The grammar of words: An introduction to linguistic morphology_.  Oxford University Press.
  * Chi et al. (2020) Ethan A Chi, John Hewitt, and Christopher D Manning. 2020.  Finding universal grammatical relations in multilingual bert.  _arXiv preprint arXiv:2005.04511_.
  * Dadegan Research Group (2012) Dadegan Research Group. 2012.  Persian dependency treebank, annotation manual and user guide.  _Supreme Council of Information and Communication Technology (SCICT), Tehran, Iran_.
  * Feely et al. (2014) Weston Feely, Mehdi Manshadi, Robert Frederking, and Lori Levin. 2014.  The CMU METAL Farsi NLP approach.  In _Proceedings of the Ninth International Conference on Language Resources and Evaluation (LREC’14)_ , pages 4052–4055, Reykjavik, Iceland. European Language Resources Association (ELRA).
  * Ghayoomi and Kuhn (2014) Masood Ghayoomi and Jonas Kuhn. 2014.  Converting an HPSG-based treebank into its parallel dependency-based treebank.  In _LREC_ , pages 802–809.
  * Grave et al. (2018) Edouard Grave, Piotr Bojanowski, Prakhar Gupta, Armand Joulin, and Tomas Mikolov. 2018.  Learning word vectors for 157 languages.  In _Proceedings of the Eleventh International Conference on Language Resources and Evaluation (LREC-2018)_ , Miyazaki, Japan. European Languages Resources Association (ELRA).
  * Hewitt and Manning (2019) John Hewitt and Christopher D. Manning. 2019.  A structural probe for finding syntax in word representations.  In _Proceedings of the 2019 Conference of the North American Chapter of the Association for Computational Linguistics: Human Language Technologies, Volume 1 (Long and Short Papers)_ , pages 4129–4138, Minneapolis, Minnesota. Association for Computational Linguistics.
  * Karimi-Doostan (2011) Gholamhossein Karimi-Doostan. 2011.  Separability of light verb constructions in persian.  _Studia Linguistica_ , 65(1):70–95.
  * Khallash et al. (2013) Mojtaba Khallash, Ali Hadian, and Behrouz Minaei-Bidgoli. 2013.  An empirical study on the effect of morphological and lexical features in Persian dependency parsing.  In _Proceedings of the Fourth Workshop on Statistical Parsing of Morphologically-Rich Languages_ , pages 97–107, Seattle, Washington, USA. Association for Computational Linguistics.
  * McDonald et al. (2013) Ryan McDonald, Joakim Nivre, Yvonne Quirmbach-Brundage, Yoav Goldberg, Dipanjan Das, Kuzman Ganchev, Keith Hall, Slav Petrov, Hao Zhang, Oscar Täckström, Claudia Bedini, Núria Bertomeu Castelló, and Jungmee Lee. 2013.  Universal dependency annotation for multilingual parsing.  In _Proceedings of the 51st Annual Meeting of the Association for Computational Linguistics (Volume 2: Short Papers)_ , pages 92–97, Sofia, Bulgaria. Association for Computational Linguistics.
  * McDonald et al. (2011) Ryan McDonald, Slav Petrov, and Keith Hall. 2011.  Multi-source transfer of delexicalized dependency parsers.  In _Proceedings of the conference on empirical methods in natural language processing_ , pages 62–72. Association for Computational Linguistics.
  * Mirzaei and Moloodi (2016) Azadeh Mirzaei and Amirsaeid Moloodi. 2016.  Persian proposition bank.  In _Proceedings of the Tenth International Conference on Language Resources and Evaluation (LREC 2016)_ , pages 3828–3835, Portorož, Slovenia. European Language Resources Association (ELRA).
  * Mirzaei and Safari (2018) Azadeh Mirzaei and Pegah Safari. 2018.  Persian discourse treebank and coreference corpus.  In _Proceedings of the Eleventh International Conference on Language Resources and Evaluation (LREC-2018)_ , Miyazaki, Japan. European Languages Resources Association (ELRA).
  * Moloodi and Kouhestani (2017) Amirsaeid Moloodi and Manouchehr Kouhestani. 2017.  The role of metaphor and metonymy in the semantics of persian adjectival preverbs: A cognitive linguistics approach.  _Language Art_ , 2(2):91–105.
  * Naseem et al. (2010) Tahira Naseem, Harr Chen, Regina Barzilay, and Mark Johnson. 2010.  Using universal linguistic knowledge to guide grammar induction.  In _Proceedings of the 2010 Conference on Empirical Methods in Natural Language Processing_ , pages 1234–1244, Cambridge, MA. Association for Computational Linguistics.
  * Nivre et al. (2020) Joakim Nivre, Marie-Catherine de Marneffe, Filip Ginter, Jan Hajič, Christopher D Manning, Sampo Pyysalo, Sebastian Schuster, Francis Tyers, and Daniel Zeman. 2020.  Universal dependencies v2: An evergrowing multilingual treebank collection.  _arXiv preprint arXiv:2004.10643_.
  * Nourian et al. (2015) Alireza Nourian, Mohammad Sadegh Rasooli, Mohsen Imany, and Heshaam Faili. 2015\.  On the importance of ezafe construction in Persian parsing.  In _Proceedings of the 53rd Annual Meeting of the Association for Computational Linguistics and the 7th International Joint Conference on Natural Language Processing (Volume 2: Short Papers)_ , pages 877–882, Beijing, China. Association for Computational Linguistics.
  * Pakzad and Minaei-Bidgoli (2016) Atefeh Pakzad and Behrouz Minaei-Bidgoli. 2016.  An improved joint model: Pos tagging and dependency parsing.  _Journal of AI and Data Mining_ , 4(1):1–8.
  * Pouramini and Mozayani (2007) Ahmad Pouramini and Naser Mozayani. 2007.  An annotation scheme for a persian treebank.  _Proceedings of Computational Linguistics In the Netherlands, CLIN_.
  * Rasooli et al. (2011a) Mohammad Sadegh Rasooli, Heshaam Faili, and Behrouz Minaei-Bidgoli. 2011a.  Unsupervised identification of Persian compound verbs.  In _Mexican International Conference on Artificial Intelligence_ , pages 394–406. Springer.
  * Rasooli et al. (2013) Mohammad Sadegh Rasooli, Manouchehr Kouhestani, and Amirsaeid Moloodi. 2013.  Development of a Persian syntactic dependency treebank.  In _Proceedings of the 2013 Conference of the North American Chapter of the Association for Computational Linguistics: Human Language Technologies_ , pages 306–314, Atlanta, Georgia. Association for Computational Linguistics.
  * Rasooli et al. (2011b) Mohammad Sadegh Rasooli, Amirsaeid Moloodi, Manouchehr Kouhestani, and Behrouz Minaei-Bidgoli. 2011b.  A syntactic valency lexicon for Persian verbs: The first steps towards Persian dependency treebank.  In _5th Language & Technology Conference (LTC): Human Language Technologies as a Challenge for Computer Science and Linguistics_, pages 227–231.
  * Rasooli and Tetreault (2015) Mohammad Sadegh Rasooli and Joel Tetreault. 2015.  Yara parser: A fast and accurate dependency parser.  _arXiv preprint arXiv:1503.06733_.
  * Rosa et al. (2014) Rudolf Rosa, Jan Mašek, David Mareček, Martin Popel, Daniel Zeman, and Zdeněk Žabokrtský. 2014.  HamleDT 2.0: Thirty dependency treebanks stanfordized.  In _Proceedings of the Ninth International Conference on Language Resources and Evaluation (LREC-2014)_ , pages 2334–2341, Reykjavik, Iceland. European Languages Resources Association (ELRA).
  * Seraji et al. (2016) Mojgan Seraji, Filip Ginter, and Joakim Nivre. 2016.  Universal dependencies for Persian.  In _Proceedings of the Tenth International Conference on Language Resources and Evaluation (LREC 2016)_ , pages 2361–2365, Portorož, Slovenia. European Language Resources Association (ELRA).
  * Seraji et al. (2014) Mojgan Seraji, Carina Jahani, Beáta Megyesi, and Joakim Nivre. 2014.  A Persian treebank with Stanford typed dependencies.  In _Proceedings of the Ninth International Conference on Language Resources and Evaluation (LREC-2014)_ , pages 796–801, Reykjavik, Iceland. European Languages Resources Association (ELRA).
  * Seraji et al. (2012) Mojgan Seraji, Beáta Megyesi, and Joakim Nivre. 2012.  Bootstrapping a Persian dependency treebank.  _Linguistic Issues in Language Technology_ , 7(18).
  * Silveira et al. (2014) Natalia Silveira, Timothy Dozat, Marie-Catherine de Marneffe, Samuel Bowman, Miriam Connor, John Bauer, and Chris Manning. 2014.  A gold standard dependency corpus for English.  In _Proceedings of the Ninth International Conference on Language Resources and Evaluation (LREC-2014)_ , pages 2897–2904, Reykjavik, Iceland. European Languages Resources Association (ELRA).
  * Straka and Straková (2017) Milan Straka and Jana Straková. 2017.  Tokenizing, pos tagging, lemmatizing and parsing ud 2.0 with udpipe.  In _Proceedings of the CoNLL 2017 Shared Task: Multilingual Parsing from Raw Text to Universal Dependencies_ , pages 88–99, Vancouver, Canada. Association for Computational Linguistics.
  * Täckström et al. (2013) Oscar Täckström, Dipanjan Das, Slav Petrov, Ryan McDonald, and Joakim Nivre. 2013.  Token and type constraints for cross-lingual part-of-speech tagging.  _Transactions of the Association for Computational Linguistics_ , 1:1–12.
  * Täckström et al. (2012) Oscar Täckström, Ryan McDonald, and Jakob Uszkoreit. 2012.  Cross-lingual word clusters for direct transfer of linguistic structure.  In _Proceedings of the 2012 Conference of the North American Chapter of the Association for Computational Linguistics: Human Language Technologies_ , pages 477–487, Montréal, Canada. Association for Computational Linguistics.
  * Taher et al. (2020) Ehsan Taher, Seyed Abbas Hoseini, and Mehrnoush Shamsfard. 2020.  Beheshti-NER: Persian named entity recognition using BERT.  _arXiv preprint arXiv:2003.08875_.
  * Tenney et al. (2019) Ian Tenney, Dipanjan Das, and Ellie Pavlick. 2019.  BERT rediscovers the classical NLP pipeline.  In _Proceedings of the 57th Annual Meeting of the Association for Computational Linguistics_ , pages 4593–4601, Florence, Italy. Association for Computational Linguistics.
  * Zeman et al. (2018) Daniel Zeman, Jan Hajič, Martin Popel, Martin Potthast, Milan Straka, Filip Ginter, Joakim Nivre, and Slav Petrov. 2018.  CoNLL 2018 shared task: Multilingual parsing from raw text to universal dependencies.  In _Proceedings of the CoNLL 2018 Shared Task: Multilingual Parsing from Raw Text to Universal Dependencies_ , pages 1–21, Brussels, Belgium. Association for Computational Linguistics.
  * Zeman and Resnik (2008) Daniel Zeman and Philip Resnik. 2008.  Cross-language parser adaptation between related languages.  In _Proceedings of the IJCNLP-08 Workshop on NLP for Less Privileged Languages_.
  * Zhang et al. (2019) Meishan Zhang, Yue Zhang, and Guohong Fu. 2019.  Cross-lingual dependency parsing using code-mixed TreeBank.  In _Proceedings of the 2019 Conference on Empirical Methods in Natural Language Processing and the 9th International Joint Conference on Natural Language Processing (EMNLP-IJCNLP)_ , pages 997–1006, Hong Kong, China. Association for Computational Linguistics.


# Neural Network-Based IRS Assisted NLoS DoA Estimation

Yasin Azhdari and Mahmoud Farhang

###### Abstract

This paper presents a learning-based approach for Direction of Arrival (DoA)
estimation using a Intelligent Reflecting Surface (IRS) or Reconfigurable
Intelligent Surface (RIS) under a Non-Line-of-Sight (NLoS) scenario. The key
innovation is the employment of a novel Neural Network (NN) layer, referred to
as the NN-based RIS layer, within a generic network structure.

The NN-based RIS layer is designed to learn the optimal RIS phase shifts that
are tailored for the DoA estimation task. To achieve this, the pre-processed
real-valued observations are fed into the RIS layer, which has a specialized
structure. Unlike regular neural network layers, the weights of the NN-based
RIS layer are constrained to be sinusoidal functions, with the phase arguments
being the tunable parameters during the training process. This allows the
layer to emulate the functionality of an RIS.

Accordingly, the standard feed-forward and back-propagation procedures are
modified to accommodate the unique structure of the NN-based RIS layer.
Numerical simulations demonstrate that the proposed machine learning-based
approach outperforms conventional non-learning-based methods for DoA
estimation under most practical SNR ranges in an RIS-assisted scheme, also it
shows a better tracking capability.

###### Index Terms:

Direction of arrival (DoA) estimation, Intelligent Reflecting Surfaces (IRSs),
Neural Network (NN).

##  I Introduction

Estimation of parameters corresponding to exponential signals contaminated in
noise is a crucial problem in various signal processing applications,
including array signal processing[1]. Array signal processing aims to estimate
parameters by exploiting both temporal and spatial information. Estimating the
angular position of some sources by a set of sensors forming an array, known
as Direction of Arrival (DoA) Estimation, is among main problems in the field
of array signal processing. DoA estimation concerns determination of the angle
of arrival of signals, in electromagnetic or acoustic wave forms, impinging on
an array of antennas. DoA estimation has a variety of applications in wireless
communications, sonar, radar, navigation and so on [2, 3, 4, 5, 6, 7].

Intelligent Reflecting Surfaces (IRSs) or Reconfigurable Intelligent Surfaces
(RISs) have been applied in numerous fields in recent years, including radars
and wireless communications. The main characteristic of the IRSs is that they
have adjustable phase, amplitude, frequency and polarization, i.e., tunable
electromagnetic response [8]. Non-line of sight paths are constructed by
proper exploit of IRS (with adjustment of its reflection coefficients) in dead
zones, where line of sight links are not established [9].

In order to enhance the performance of a detection system, an IRS-assisted
system is proposed which adjusts the IRS phases to improve the received
signals’ SNR [10]. Moreover, the phase of the IRS is adjusted to maximize the
signal-to-noise ratio (SNR) in the direction corresponding to the user
equipment for localization in wireless communication under near-field
propagation regime [11]. In addition, the IRS-assisted scenario is employed
for Direction of Arrival (DoA) estimation problem, too. For example, in [12],
a coprime linear array (CLA) is implemented by controlling IRS units, also a
DoA estimation method is presented corresponding to the IRS-based CLA. In
[13], a cost-effective direction-finding system using an unmanned aerial
vehicle (UAV) swarm is presented. This system includes a central full-
functional receiving unit that performs DoA estimation by solely receiving the
signals reflected by the IRS. Furthermore, [14] explores a method where the
array receiver collects both the reflected signal from the IRS and the direct-
path signal. By designing the phase of the IRS, the reflected path can be
leveraged to enhance the accuracy of DoA estimation. Additionally, [15]
formulates and addresses the DoA estimation problem in the presence of
wireless communication interference using an IRS. An atomic norm-based
approach is then proposed for joint DoA estimation and interference removal,
where the optimization problem is also solved to design the measurement matrix
(the matrix containing IRS characteristics) for interference mitigation.

An important subject under IRS-related topics is how to adjust and tune its
responce (specifically the phase response), i.e., measurement matrix, in order
to attain the desired performance. In related works phase design is
implemented to maximize the SNR towards the desired user equipment, to
maximize the coherence of the signals from direct and reflected paths, to
minimize the Cramer-Rao Lower Bound (CRLB) of the problem, or to remove the
interference [11, 14, 16, 15]. Though a reasonable way is to use Artificial
Neural Networks (ANNs) to optimize the RIS phases, as far as we know, no
particular task-based phase design using a dedicated Neural Network (NN) layer
as the RIS layer is proposed yet.

An Artificial Neural Network (ANN) can approximate complex functions, it has
also shown the capability to establish the mapping between some input features
(usually the raw measurements or the correlation matrix) and directions of
signals’ sources [17]. For example, Multi layer perceptron (MLP) neural
networks were used for the problem of DoA estimation [18]. They are usually
constructed from a number of hidden layers with corresponding activation
functions. For instance [17] presented a direction-of-arrival estimation
method with high-resolution for the case of a single snapshot. The literature
[19] proposed to use real and imaginary components of the covariance matrix as
two seperate channels for the input of the Convolutional Neural Network (CNN)
in order to estimate direction of arrivals, and therefore the DoA estimation
problem was transformed into an image classification problem. In addition,
[20] proposed a combination of a linear classifier network (LCN) and a
convolutional neural network, where the angular range was divided into smaller
sectors by the former, then each sector included one CNN for the DoA
estimation purpose. Moreover, [21] proposed a scheme combining regression and
classification in order to simultaneously estimate DoAs and the number of
sources.

So the main contribution of this paper is introduction of a novel layer for
ANNs, which acts as the RIS (NN-based RIS layer) and providing its
corresponding mathematics in order to learn the optimal phase design with
respect to the task of DoA estimation, which has not been done yet to our best
knowledge. Also we will derive the Cramer-Rao Lower Bound (CRLB) and also
Maximum Likelihood (ML) estimation corresponding to the problem model and
formulation respecting previous works and Compare the performance of the
proposed learning-based method with the ML estimate and also non learning
methods modified for our specific problem model and formulation. In this
research, Uniform Linear Arrays (ULAs) are considered. But the method can be
easily extended to other types of arrays.

The remainder of this paper is structured as follows: The problem model and
formulation is described in Section II. Section III reviews well-known non
learning-based RIS assisted DoA estimation techniques and their corresponding
algorithms. Section IV presents the novel NN-based RIS layer and its
corresponding mathematics, employed in generic neural network structures for
the task of DoA estimation; numerical simulations and results are provided in
Section V. Finally, conclusions are given in Section VI.

Notations: $j=\sqrt{-1}$ represents the imaginary unit. Boldface lower case
letters refer to vectors and boldface upper case letters denote matrices. For
a square matrix $\mathbf{S}$, $\mathbf{S}^{-1}$ and tr($\mathbf{S}$) denote
its inverse and trace, respectively. For an arbitrary matrix $\mathbf{M}$,
$\mathbf{M}^{T}$, $\mathbf{M}^{H}$, $\mathbf{M}^{*}$, and rank($\mathbf{M}$)
are its transpose, conjugate transpose, conjugate, and rank, respectively.
Also $\mathbf{M}^{+}$ denotes its Moore-Penrose or pseudo inverse defined as
$(\mathbf{M}^{T}\mathbf{M})^{-1}\mathbf{M}^{T}$. The matrix $\mathbf{I}_{m}$
represents the identity matrix with dimension $m\times m$. We use
$CN(\mathbf{m};\mathbf{C})$ to show the distribution of a complex Gaussian
random vector with mean of vector $\mathbf{m}$ and covariance of matrix
$\mathbf{C}$. The real and imaginary parts of a complex number are denoted by
Re{.} and Im{.}, respectively. E{.} stands for the statistical expectation.
$|.|$ represents the magnitude of a complex number and $\|.\|$ stands for the
Euclidean norm. $diag(a_{1},...,a_{N})$ represents the diagonalization
operator which results in a diagonal matrix with diagonal elements of
$a_{1},...,a_{N}$. $\odot$ denotes element-wise product between two matrices
or vectors of the same dimension.

##  II Problem Model

###  II-A Formulation

We consider an IRS-assisted scenario as shown in Fig. 1. According to figure
1, it is assumed that there is no line of sight path between the source and
the antenna array. Thus the only path available is through IRS reflection. The
antenna array and the IRS are both considered to be ULAs. The antenna array is
placed on the z axis of the coordinate system and the IRS is placed on the x
axis. The intersection of the corresponding axes is assumed to be the origin
of the Cartesian coordinate system. The target is assumed to be in the far-
field region with respect to the IRS and antennas, whereas the IRS is placed
in the near-field region with respect to the antennas. Without loss of
generality, this paper assumes one target at the non line of sight far-field
region of the array, i.e., the most powerful one. Narrow-band waveform is
assumed to prevent the effect of signal propagation on the IRS and the array.
Also, unit gains and isotropic power receiving patterns are considered for
both array antennas and IRS unit cells. As optimizing those factors are out of
this paper’s scope.

Figure 1: Problem Model

$M^{A}$ and $M^{R}$ denote the number of ULA elements and IRS unit cells, with
the inter element spacing of $d^{A}$ and $d^{R}$, respectively. The angle
between the wavefront and the orthogonal line with respect to the RIS axis is
denoted as $\theta$. ${p_{m}}^{A}=({x_{m}}^{A},{z_{m}}^{A})$ is the position
of the mth array element. ${p_{n}}^{R}=({x_{n}}^{R},{z_{n}}^{R})$ and
$\gamma_{n}=|\gamma_{n}|e^{j\phi_{n}}$ denote the position and the adjustable
reflection coefficient of the nth unit cell, respectively. $|\gamma_{n}|$s are
assumed to be 1 throughout this paper.
${r_{m}}^{AT}=(m-1)d^{A}\sin(\theta_{m,n}^{AR})$ denotes the path difference
between the mth array element and the reference phase element.
${r_{n}}^{RT}=(n-1)d^{R}\sin(\theta)$ represents the path difference between
the reference phase cell and the nth IRS unit cell.
${r_{m,n}}^{AR}=\|{p_{m}}^{A}-{p_{n}}^{R}\|_{2}$ and ${\theta_{m,n}}^{AR}$
denotes the distance between the nth IRS unit cell and the mth array element
and the corresponding angle shown in Fig. 1, respectively.

Considering $\mathbf{\Phi}$ as the IRS phase shifts vector,
$\mathbf{\omega}=e^{j\mathbf{\Phi}}=[e^{j\phi_{1}},e^{j\phi_{2}},...,e^{j\phi_{M^{R}}}]$
denotes the corresponding reflective steering vector. Then $\Omega$ is defined
as the $M^{R}\times M^{R}$ matrix containing the IRS phases as follows:

$\mathbf{\Omega}=\text{diag}(\mathbf{\omega})=\text{diag}([e^{j\phi_{1}},e^{j\phi_{2}},\ldots,e^{j\phi_{M^{R}}}])$
(1)

Considering a far-field deterministic source impinging, the recieved signal
from the reflection path can be modeled as [11, 14]:

$\mathbf{y_{r}}=(\mathbf{H}^{AR}\odot\mathbf{A}^{AR})\mathbf{\Omega}\mathbf{a}^{RT}s=\mathbf{a}_{r}s$
(2)

where $\mathbf{A}^{AR}$ and $\mathbf{a}^{RT}$ are the array-IRS steering
matrix and IRS-target steering vector, respectively.

$\mathbf{A}^{AR}=\left[\begin{array}[]{lr}e^{-j2\pi
r_{1,1}^{AR}/\lambda}...&e^{-j2\pi r_{1,M^{R}}^{AR}/\lambda}\\\ ...&...\\\
e^{-j2\pi r_{M^{A},1}^{AR}/\lambda}...&e^{-j2\pi
r_{M^{A},M^{R}}^{AR}/\lambda}\end{array}\right]$ (3)

$\mathbf{a}^{RT}=[e^{j2\pi r_{1}^{RT}/\lambda},...,e^{j2\pi
r_{M^{R}}^{RT}/\lambda}]^{T}$ (4)

Also, $\mathbf{H}^{AR}$ is the matrix representing the amplitude of the
received signal from the near-field reflection paths [22]

${H}_{m,n}^{AR}=(\frac{P_{r}\lambda^{2}}{(4\pi{r_{m,n}}^{AR})^{2}})^{\frac{1}{2}}$
(5)

where $P_{r}$ denotes the power of the signal approaching the IRS and the
array receiver, which is assumed to be equal due to the far-field assumption.
This $P_{r}$ value is equal to $P_{t}-\text{Path loss}_{\text{target-IRS}}$
(in dB), where $P_{t}$ is the transmitted signal power, and the Path loss is
defined as follows:

$L(d)=K_{0}(\frac{d}{d_{0}})^{-\alpha_{0}}$ (6)

where $d$ represents the length of the transmission link, $K_{0}$ denotes the
path loss at the length = 1 m, $\alpha_{0}$ is the corresponding path loss
exponent.

Thus the recieved signal at array can be written as:

$\mathbf{y}=\mathbf{y_{r}}+\mathbf{n}=(\mathbf{H}^{AR}\odot\mathbf{A}^{AR})\mathbf{\Omega}\mathbf{a}^{RT}s+\mathbf{n}=\mathbf{a}_{r}s+\mathbf{n}$
(7)

where $\mathbf{n}$ is the additive noise distributed as
$CN(\mathbf{0},\sigma_{n}^{2}I)$.

###  II-B Cramer-Rao Lower Bound

The CRLB for our specific problem model and formulation, can be obtained in a
similar fashion to [23]. The likelihood function and its corresponding
logarithm of the received measurements at array is given by (signal and noise
variance are assumed to be known):

$\text{L}_{y}=\frac{1}{(\pi\sigma^{2})^{M^{A}}}\text{exp}(-\frac{1}{\sigma^{2}}(\mathbf{y}-\mathbf{a}_{r}(\theta)s)^{H}(\mathbf{y}-\mathbf{a}_{r}(\theta)s))$
(8)

$\text{ln}\text{L}_{y}=\text{const}-\frac{1}{\sigma^{2}}(\mathbf{y}-\mathbf{a}_{r}(\theta)s)^{H}(\mathbf{y}-\mathbf{a}_{r}(\theta)s)$
(9)

The Fisher information measure and the CRLB corresponded to DoA estimation
problem in the aforementioned RIS-assisted scheme are defined as:

$\text{CRLB}_{\theta}=\text{FIM}^{-1}(\theta)=-\text{E}\\{\frac{\partial^{2}}{\partial\theta^{2}}\text{ln}\text{L}_{y}\\}^{-1}$
(10)

Inserting equation 9 in 10 results in $\textit{CRLB}_{\theta}$ to be:

$\frac{\sigma_{n}^{2}}{2\sigma_{s}^{2}}\text{Re}\\{\frac{\partial\mathbf{a}_{r}(\theta)}{\partial\theta}^{H}(I-\mathbf{a}_{r}(\theta)(\mathbf{a}_{r}(\theta)^{H}\mathbf{a}_{r}(\theta))^{-1}\mathbf{a}_{r}(\theta)^{H})\frac{\partial\mathbf{a}_{r}(\theta)}{\partial\theta}\\}^{-1}$
(11)

where $\sigma_{s}^{2}=\text{E}\\{\|s\|^{2}\\}=\|s\|^{2}$ and

$\frac{\partial\mathbf{a}_{r}(\theta)}{\partial\theta}=(\mathbf{H}^{AR}\odot\mathbf{A}^{AR})\mathbf{\Omega}\frac{\partial\mathbf{a}^{RT}(\theta)}{\partial\theta}$
(12)

where

$\frac{\partial\mathbf{a}^{RT}(\theta)}{\partial\theta}=\frac{\partial}{\partial\theta}[e^{j2\pi
r_{1}^{RT}/\lambda},...,e^{j2\pi r_{M^{R}}^{RT}/\lambda}]^{T}$ (13)

considering that ${r_{n}}^{RT}=(n-1)d^{R}\sin(\theta)$, (13) turns in to

$[0,\ldots,\frac{j2\pi(M^{R}-1)d^{R}\cos(\theta)}{\lambda}e^{j2\pi(M^{R}-1)d^{R}\sin(\theta)/\lambda}]^{T}$
(14)

Thus the corresponding CRLB of the problem model and formulation mentioned
before, which determines a lower bound on the Mean Squared Error (MSE) of any
unbiased estimators for $\theta$, can be obtained accordingly.

##  III Non learning-based RIS assisted DoA estimation

###  III-A Phase adjustment method

Considering the fact that in the aforementioned non line of sight scenario,
maximizing the coherence between direct and reflected path signals can not be
considered, and also due to the fact that the phase adjustment based on
maximizing the received SNR at the array antenna is outperformed by CRLB
minimization criterion [16], we will consider minimizing CRLB by selecting
proper set of phases $\mathbf{\Phi}$ in our specific problem model.

$\text{min}_{\mathbf{\omega}}\;J(\mathbf{\omega})=\text{CRLB}_{\theta}(\mathbf{\omega})\;\;\text{s.t.}\;\;|\omega_{i}|=1\;\forall
i=1,...,M^{R}$ (15)

where $\text{CRLB}_{\theta}(\omega)$ is given in (11). Because the CRLB
function and the unit-modulus constraints are not convex, The optimization
problem of the IRS phase design is not convex and thus is challenging and
complicated. In order to transform this problem into a convex optimization
problem, the semi definite relaxation (SDR) and successive convex
approximation (SCA) techniques can be employed [16]. Also the method of
Riemannian manifold optimization can be applied for designing optimal phase
shifts [14]. Here we have complex circle manifold, and to solve the Riemannian
manifold-based optimization problem (15) the Riemannian Steepest Descent (SD),
Conjugate Gradient (CG), and Trust-Region (TR) methods are proposed before
[24].

Brief detail of the corresponding algorithms are presented in the Algorithm 1.
Here we will utilize the well-known TR method.

Algorithm 1 Manifold-based optimization for CRLB minimization

Initialize $\mathbf{\omega}_{0}\in\mathcal{M}$, $i=0,\epsilon$

while $|\nabla_{\mathcal{M}}J(\omega_{i})|<\epsilon$ do

1) Compute the corresponding Riemannian gradient of the cost function
$\nabla_{\mathcal{M}}J(\omega_{i})$ on the complex circle manifold
$\mathcal{M}$;

2) Compute the corresponding Riemannian hessian of the cost function
$\nabla_{\mathcal{M}}J(\omega_{i})$ on the complex circle manifold
$\mathcal{M}$ (in case needed, i.e., for Riemannian TR solver);

3) update the solution and get $\omega_{i+1}$ based on the Riemannian SD, CG,
or TR algorithms [24];

4) $i=i+1$

end while

###  III-B Doa estimation methods

####  III-B1 Maximum Likelihood Doa Estimation

The Maximum Likelihood (ML) DoA estimator for the deterministic signal is
implemented as follows and requires a grid search [25]:

$\hat{\theta}=\textit{argmax}_{\theta}\;\lambda_{max}[(\mathbf{a}_{r}^{H}\mathbf{a}_{r})^{-\frac{1}{2}}(\mathbf{a}_{r}^{H}yy^{H}\mathbf{a}_{r})(\mathbf{a}_{r}^{H}\mathbf{a}_{r})^{-\frac{1}{2}}]$
(16)

where $\lambda_{max}[\mathbf{X}]$ denotes the maximum eigenvalue of the matrix
$\mathbf{X}$.

####  III-B2 MUSIC

The MUSIC (Multiple Signal Classification) algorithm [26] is a widely-used
subspace-based direction-of-arrival (DoA) estimation technique. It operates by
separating the orthogonal noise and signal subspaces, which are derived from
the eigenvectors of the covariance matrix
($\mathbf{R}=E\\{\mathbf{y}\mathbf{y}^{H}\\}$). Specifically, if
$\mathbf{v}_{1}$ to $\mathbf{v}_{M}$ are the eigenvectors of $\mathbf{R}$,
sorted in descending order of their corresponding eigenvalues, and $K$ is the
number of sources (which here is assumed to be one), then $\mathbf{v}_{1}$ to
$\mathbf{v}_{K}$ and their linear combinations define the signal subspace,
while $\mathbf{v}_{K+1}$ to $\mathbf{v}_{M}$ and their linear combinations
define the noise subspace. Exploiting the orthogonality between the signal and
noise subspaces, the MUSIC output power is then defined as:

$P(\theta)=\frac{1}{\sum_{i=K+1}^{M}|\mathbf{a}_{r}^{H}\mathbf{v}_{i}|^{2}}=\frac{1}{\mathbf{a}_{r}^{H}\mathbf{V}_{n}\mathbf{V}_{n}^{H}\mathbf{a}_{r}}$
(17)

where $\mathbf{v}_{K+1}$ to $\mathbf{v}_{M}$ form the matrix $\mathbf{V}_{n}$.

####  III-B3 Compressed Sensing (CS)-based Doa Estimation

The sparse reconstruction problem for the multiple measurement vector case is
expressed as follows [27]:

$\text{min}_{\mathbf{x}}\;\frac{1}{2}\|\mathbf{y-A_{all}x}\|_{2,2}^{2}+\tau\|\mathbf{x}\|_{2,1}$
(18)

where the over completed dictionary $\mathbf{A_{all}}$ of size $M^{A}\times
M^{D}$ is given as
$\mathbf{A_{all}}=[\mathbf{a}(\theta_{1}),\mathbf{a}(\theta_{2}),...,\mathbf{a}(\theta_{M^{D}})]$,
the sparse spatial spectrum is denoted as $\mathbf{x}$, and the non negative
regularization parameter is represented as $\tau$ which can be determined by
[27]:

$\tau=\sqrt{\frac{M^{D}L\sigma_{n}^{2}}{K}}$ (19)

where $L$ , $K$ denote the number of snapshots and the number of sources
(which here is assumed to be 1), respectively. The convex problem (18) can be
tackled by interior-point algorithm-based packages, i.e., CVX [28].

##  IV Learning-based IRS assisted DoA estimation

As mentioned earlier, in order to adjust IRS units phases, several non
learning-based optimization methods have been proposed. In this paper we have
implemented the IRS using a dedicated novel neural network layer, which learns
to adjust the phases considering the task of DoA estimation, and presented its
corresponding mathematics. Since the overall network will use raw measurements
observed at the IRS as inputs, and because inputs of the network should be
real-valued, first we pre-process them by seperating the real and imaginary
parts of each unit’s measurement and also normalizing them. Then the pre-
processed observations are fed to the IRS layer. But the IRS layer can not be
the same as other ordinary layers with tunable constant weights and its
weights should be sinusoidals, where the arguments or phases of these
sinusoidals should be tuned during the training procedure. So we need to
modify the feed forward and backward propagation procedures’ mathematics
accordingly. Moreover, the IRS layer will utilize the linear activation
function. In the following we will clarify the corresponding process.

###  IV-A Implementation of the IRS layer

A fully connected (FC) layer operates based on forward propagation and
backward propagation procedures. Considering a FC layer with input vector
$\mathbf{x}\in R^{1\times m}$, weights matrix $\mathbf{W}\in R^{n\times m}$,
bias vector $\mathbf{b}\in R^{1\times n}$, and output vector $\mathbf{y}\in
R^{1\times n}$, forward propagation (FP) formulation is as following:

$\textit{FP}:\mathbf{y}=\mathbf{xW}+\mathbf{b}$ (20)

Then the output of one layer becomes the input of the next according to
equation (20). Then, considering the result and the desired output, an error
$E$ is obtained per each output neuron. The goal is to minimize that error by
changing the parameters in the network, i.e., backward propagation based on
the gradient decsent algorithm. In order to do such, assuming that we have a
matrix including the derivative of the error respecting the layer’s output
($\frac{\partial E}{\partial\mathbf{y}}$), it is required to obtain the
derivative of the error respecting parameters ($\frac{\partial
E}{\partial\mathbf{W}},\frac{\partial E}{\partial\mathbf{b}}$) and also the
input ($\frac{\partial E}{\partial\mathbf{x}}$) using the chain rule.

Next, we will implement forward propagation and backward propagation
procedures dedicated for our specific desired IRS layer. First we define the
weights matrix to be as follows:

$\mathbf{W}=\left[\begin{array}[]{lr}\gamma_{1,1}...&\gamma_{1,m}\\\
...&...\\\ \gamma_{n,1}...&\gamma_{n,m}\end{array}\right]\in R^{n\times m}$
(21)

where $\gamma_{i,j}=\sin(\theta_{i,j})\;\forall i=1,...,n,j=1,...,m$.

The logic behind defining the weights matrix as (21) comes from the fact that
we want the IRS layer to multiply exponentials with optimum phase arguments in
observations, resulting in an optimal linear combination of observations with
exponential weights. Considering that we have seperated real and imaginary
parts of observations in the pre-processing phase, here the exponentials turn
into sinusoidals, where optimum phases (arguments) are learned during the
training phase. With weights matrix defined, next we will derive the
aforementioned derivatives. The $\frac{\partial E}{\partial\mathbf{W}}$ matrix
is defined as follows:

$\frac{\partial E}{\partial\mathbf{W}}=\left[\begin{array}[]{lr}\frac{\partial
E}{\partial W_{1,1}}...&\frac{\partial E}{\partial W_{1,m}}\\\ ...&...\\\
\frac{\partial E}{\partial W_{n,1}}...&\frac{\partial E}{\partial
W_{n,m}}\end{array}\right]\in R^{n\times m}$ (22)

But since we want to adjust and tune $\theta_{i,j}$s and not ${W}_{i,j}$s, we
will use chain rule once more as follows:

$\frac{\partial E}{\partial\theta_{i,j}}=\sum_{j=1}^{m}\frac{\partial
E}{\partial y_{j}}.\frac{\partial
y_{j}}{\partial\gamma_{i,j}}.\frac{\partial\gamma_{i,j}}{\partial\theta_{i,j}}$
(23)

Since $y_{j}=b_{j}+\sum_{i=1}^{n}x_{i}\gamma_{i,j}$, $\frac{\partial
y_{j}}{\partial\gamma_{i,j}}$ equals zero except for one term, thus:

$\frac{\partial E}{\partial\theta_{i,j}}=\frac{\partial E}{\partial
y_{j}}.x_{i}.\cos(\theta_{i,j})$ (24)

So $\frac{\partial E}{\partial\theta}$ can be formulated as follows:

$\begin{split}\frac{\partial
E}{\partial\mathbf{\theta}}&=\text{diag}(\mathbf{x}).\left[\begin{array}[]{lr}\cos(\theta_{1,1})...&\cos(\theta_{1,m})\\\
...&...\\\
\cos(\theta_{n,1})...&\cos(\theta_{n,m})\end{array}\right].\text{diag}(\frac{\partial
E}{\partial\mathbf{y}})\\\ &=\mathbf{x}^{D}.\mathbf{C}.(\frac{\partial
E}{\partial\mathbf{y}})^{D}\end{split}$ (25)

Next, we will consider the $\frac{\partial E}{\partial\mathbf{x}}$. Since:

$\frac{\partial E}{\partial x_{i}}=\sum_{j=1}^{m}\frac{\partial E}{\partial
y_{j}}.\frac{\partial y_{j}}{\partial x_{i}}=\sum_{j=1}^{m}\frac{\partial
E}{\partial y_{j}}.\gamma_{i,j}$ (26)

Thus:

$\frac{\partial E}{\partial\mathbf{x}}=\frac{\partial
E}{\partial\mathbf{y}}.\mathbf{W}^{T}$ (27)

Additionaly, here we haven’t considered bias for the IRS layer.

Finally, since we want the output to be a linear combination of inputs with
sinusoidal weights, the corresponding activation for the IRS layer is
considered to be linear.

The gradient descent algorithm with proper learning factor is applied to
update the weight parameters $\theta_{i,j}$ in order to minimize the loss
function, which considering the nature of DoA regression problem has chosen to
be the Mean Squared Error (MSE), defined as follows:

$E=\frac{1}{n}\sum_{i=1}^{n}({y_{i}}^{*}-y_{i})^{2}$ (28)

with the following derivative respecting the layer’s output:

$\frac{\partial E}{\partial\mathbf{y}}=\frac{2}{n}(\mathbf{y}-\mathbf{y}^{*})$
(29)

###  IV-B Implementation of the network

The input of the network is the observations of the IRS unit cells. Thus in
order to formulate the network’s input we consider the target-IRS link as a
direct path and formulate IRS observations as the following:

$\mathbf{x}=\mathbf{a}^{RT}s+\mathbf{n}$ (30)

where $\mathbf{x}$ denotes the IRS unit cells’ observations. Also,
$\mathbf{a}^{RT}$ and $\mathbf{n}$ represent the IRS-target steering vector
and the additive noise as discussed before. The loss corresponding to the
target-IRS link path is also taken into the account. Next, we have the novel
NN-based IRS layer with tunable phases as discussed before. After that we will
utilize two different network architectures, including a fully connected and
also a convolutional network. Moreover, a residual connection will also be
considered in order to maintain an additive direct link from the input to the
DoA regressor network beginning from just after the IRS layer.

Assuming a RIS with 10 unit cells and under 10 snapshots, scheme of the
overall fully connected and also convolutional networks are illustrated in
figures 2, 3, respectively. Which represent both the input size and the output
size of each layer and also corresponding activation functions.

Figure 2: The overall fully connected network Figure 3: The overall
convolutional network

Both networks include a pre-processing block, which turns each IRS unit’s
observation into a 2-dimensional vector containing its real and imaginary
parts and normalizes them. Then we have the novel IRS block described in the
previous subsection, calculating a linear combination of the inputs with
sinusoidal weights, whose phases are to be learned considering the task of DoA
estimation. Next we got a supervised DoA regressor network implemented under
two different architectures. Thus the overall network can be represented with
the following transform $\text{N}(.)$:

$\hat{\theta}=\text{N}(\mathbf{z})=H_{9}(H_{8}(...(H_{1}(R(P(\mathbf{z})))))):R^{10\times
10}\rightarrow R$ (31)

where $P(.)$ represents the pre-processing performed. Moreover,
$R(\mathbf{u})$ represents the novel IRS layer where $\mathbf{W}_{r}$
represents its specific weights matrix defined in (21). Additionaly
$H_{i}(\mathbf{u})$s are representing hidden layers of the networks where
$\mathbf{W}_{i}$s represent their weights matrix. Finally the output of the
network is the estimated DoA approximated by the neural network structure.

The overall network, including both the IRS block and the DoA estimator block
will be tuned and adjusted simultaneously during the training phase concerning
the ultimate task of DoA estimation. The details and parameters of the network
and also the generation procedure of the corresponding datasets are presented
in the next section, where the performance of the proposed scheme is compared
with other IRS-based DoA estimation methods for the same problem model through
simulations.

##  V Numerical Simulations and Results

###  V-A Setup

The number of array elements and the number of IRS unit cells are assumed to
be $M^{A}=M^{R}=10$. Problem model with simulation presets is illustrated in
Fig. 4. Source signal is assumed to be an narrow-band deterministic
exponential with frequency of 1 GHz. Since the far-field range
$>\frac{2D^{2}}{\lambda}$ (where $D$ denotes the largest antenna dimension)
and considering that for both the IRS and the array with inter-element spacing
of $\lambda/2$, $D=(M^{R}-1)\lambda/2=(M^{A}-1)\lambda/2=2\lambda=0.6m$, the
far-field range in our considered setup will be any distance greater than 2.4
meters. Based on the above far-field range threshold, we have considered the
distances mentioned in Fig. 4. $K_{0}$ =-15 dB is the path loss at the
reference distance $d_{0}$ = 1 m, and the path-loss exponent $\alpha_{0}$ is
set to be 2.5. the number of snapshots is 1 and the noise variance is assumed
to be $\sigma_{n}^{2}=1$. The SNR is defined as
$10\log_{10}\frac{P_{t}}{\sigma_{n}^{2}}=10\log_{10}(\frac{|s|^{2}}{\sigma_{n}^{2}})$.
The network’s input and output dimensions, activations, and also the number of
parameters per each layer are gathered in Table 1. The total number of
parameters are 1,708,501 for the fully connected network and 2,955,453 for the
convolutional one.

(0,0)1 m1 m10 m$\theta$

Figure 4: Problem Model with Simulation Presets TABLE I: Neural Network details Layer | Input dim. (F.C.) | Output dim. (F.C.) | Input dim. (Conv.) | Output dim. (Conv.) | Activation
---|---|---|---|---|---
$P(.)$ | (10,10) | (10,20) | (10,10) | (10,10) | -
$R(.)$ | (10,20) | (10,20) | (10,20) | (10,20) | linear
$H_{1}(.)$ | 200 | 200 | (10,20,1) | (8,18,16) | tanh
$H_{2}(.)$ | 200 | 250 | (8,18,16) | (6,16,32) | tanh
$H_{3}(.)$ | 250 | 500 | (6,16,32) | (4,14,64) | tanh
$H_{4}(.)$ | 500 | 750 | (4,14,64) | (2,12,128) | tanh
$H_{5}(.)$ | 750 | 750 | 3072 | 750 | leaky relu
$H_{6}(.)$ | 750 | 500 | 750 | 500 | leaky relu
$H_{7}(.)$ | 500 | 250 | 500 | 250 | leaky relu
$H_{8}(.)$ | 250 | 50 | 250 | 50 | leaky relu
$H_{9}(.)$ | 50 | 1 | 50 | 1 | linear

The output values (DoAs) are normalized to the largest Field of View (FoV)
value during training and test phases. Moreover, as described before, the loss
function is assumed to be the RMSE.

Initialization for both manifold-based optimization algorithm’s phases and
also the RIS layer’s sinusoidal weights phases is done based on a uniform
distribution over the set $[-\pi,\pi]$.

###  V-B DataSets

Our training set includes $N_{train}=10,000$ examples, i.e. an input-DoA pair
$(\mathbf{x},\theta)$. the number of DoA sources is assumed to be one and 10
snapshots are considered. Source DoAs are distributed uniformly over the FoV,
which contains from $-45^{\circ}$ to $45^{\circ}$. The SNR values are drawn
independently according to a uniform distribution over the set of
$\\{-10,-5,0,5,10,15,20\\}$ dBs.

The test set contains $N_{test}=100$ examples per each SNR, then by averaging
the evaluated performance over $100$ independent test realizations per each
SNR value, the results are obtained.

###  V-C Experiments

Assuming a source at $\theta=12.5^{\circ}$, the root mean squared error (RMSE)
for proposed learning based methods (Fully connected and Convolutional
networks) and also the non-learning-based ones (ML, MUSIC, and Compressed
Sensing (CS)) are plotted versus SNR through figures 6, 7, 8, and 9. The Root
Mean Squared Error (RMSE) is defined as:

$\mbox{RMSE}=\sqrt{\frac{1}{C}\sum_{c=1}^{C}(\hat{\theta}_{c}-\theta)^{2}}$
(32)

where $C$ denotes number of MonteCarlo simulation runs (considered 100 times),
$\theta$ denotes angular location of the source and $\hat{\theta}_{c}$
represents its estimation resulted during $c$-th run (in degrees).

Figure 5: RMSE versus SNR curves - Fully connected networks Figure 6: RMSE
versus SNR curves - Convolutional networks Figure 7: RMSE versus SNR curves -
Fully connected vs. Convolutional Figure 8: RMSE versus SNR curves - Learning
vs. Non-learning

As illustrated by the line graphs, for the case of convolutional networks it
is better to utilize the residual connection, while for the case of fully
connected one it isn’t. Moreover, for low SNR values (below 0) the
convolutional network performs slightly better, while for higher SNR values
the fully connected network outperforms it. The proposed neural network
(NN)-based RIS-assisted Direction of Arrival (DoA) estimator outperforms the
Compressed Sensing (CS) technique based on manifold optimization for almost
all SNRs. Also it outperforms the MUSIC approach based on manifold
optimization under low SNR values, which is the case in most practical
scenarios.

Next, considering a source with constant angular velocity, moving from
$-45^{\circ}$ to $45^{\circ}$, the averaged spatial spectrums over 100
independent runs are depicted in Fig. 9 which confirms that the proposed
method has lower fluctuations while tracking the source and performs better.

Figure 9: Average Predicted angle versus Actual angle

##  VI Conclusion

This paper introduced a novel approach that utilizes a dedicated fully
connected neural network layer to act as the Reconfigurable Intelligent
Surface (RIS) (referred to as the NN-based RIS layer). This novel NN-based RIS
layer is designed to learn the optimal phase configuration for the purpose of
Direction of Arrival (DoA) estimation within a generic multilayer perceptron
(MLP) structure.

The proposed method performs the RIS phase adjustment and DoA estimation tasks
simultaneously, without requiring any additional optimization steps. Numerical
simulation results demonstrate the superior performance of this approach
compared to other non-learning-based DoA estimation techniques in an RIS-
assisted scheme under most practical SNR values and also under a tracking
scenario.

## References

  * [1] H. Ye and D. DeGroat, “Maximum likelihood doa estimation and asymptotic cramér-rao bounds for additive unknown colored noise,” _IEEE Transactions on Signal Processing_ , vol. 43, no. 4, pp. 938–949, 1995.
  * [2] H. Krim and M. Viberg, “Two decades of array signal processing research: the parametric approach,” _IEEE signal processing magazine_ , vol. 13, no. 4, pp. 67–94, 1996.
  * [3] H. L. Van Trees, _Optimum array processing: Part IV of detection,estimation, and modulation theory_. Wiley-Interscience, 2008.
  * [4] P.-J. Chung, M. Viberg, and J. Yu, “Doa estimation methods and algorithms,” in _Academic Press Library in Signal Processing_. Elsevier, 2014, vol. 3, pp. 599–650.
  * [5] R. O. Nielsen, _Sonar signal processing_. Artech House, Inc., 1991.
  * [6] Z. Chen, G. Gokeda, and Y. Yu, _Introduction to Direction-of-arrival Estimation_. Artech House, 2010.
  * [7] T. E. Tuncer and B. Friedlander, _Classical and modern direction-of-arrival estimation_. Academic Press, 2009.
  * [8] C. Huang, A. Zappone, G. C. Alexandropoulos, M. Debbah, and C. Yuen, “Reconfigurable intelligent surfaces for energy efficiency in wireless communication,” _IEEE transactions on wireless communications_ , vol. 18, no. 8, pp. 4157–4170, 2019.
  * [9] D. Dardari, N. Decarli, A. Guerra, and F. Guidi, “Los/nlos near-field localization with a large reconfigurable intelligent surface,” _IEEE Transactions on Wireless Communications_ , vol. 21, no. 6, pp. 4282–4294, 2021\.
  * [10] S. Buzzi, E. Grossi, M. Lops, and L. Venturino, “Radar target detection aided by reconfigurable intelligent surfaces,” _IEEE Signal Processing Letters_ , vol. 28, pp. 1315–1319, 2021.
  * [11] A. Elzanaty, A. Guerra, F. Guidi, and M.-S. Alouini, “Reconfigurable intelligent surfaces for localization: Position and orientation error bounds,” _IEEE Transactions on Signal Processing_ , vol. 69, pp. 5386–5402, 2021.
  * [12] T. Lan, K. Huang, L. Jin, X. Xu, X. Sun, and Z. Zhong, “Doa estimation algorithm for reconfigurable intelligent surface co-prime linear array based on multiple signal classification approach,” _Information_ , vol. 13, no. 2, p. 72, 2022.
  * [13] P. Chen, Z. Chen, B. Zheng, and X. Wang, “Efficient doa estimation method for reconfigurable intelligent surfaces aided uav swarm,” _IEEE Transactions on Signal Processing_ , vol. 70, pp. 743–755, 2022.
  * [14] H. Chen, Y. Bai, Q. Wang, H. Chen, L. Tang, and P. Han, “Doa estimation assisted by reconfigurable intelligent surfaces,” _IEEE Sensors Journal_ , 2023.
  * [15] Z. Chen, P. Chen, Z. Guo, Y. Zhang, and X. Wang, “A ris-based vehicle doa estimation method with integrated sensing and communication system,” _IEEE Transactions on Intelligent Transportation Systems_ , 2023.
  * [16] X. Song, J. Xu, F. Liu, T. X. Han, and Y. C. Eldar, “Intelligent reflecting surface enabled sensing: Cramér-rao bound optimization,” _IEEE Transactions on Signal Processing_ , 2023.
  * [17] J. Fuchs, R. Weigel, and M. Gardill, “Single-snapshot direction-of-arrival estimation of multiple targets using a multi-layer perceptron,” in _2019 IEEE MTT-S International Conference on Microwaves for Intelligent Mobility (ICMIM)_. IEEE, 2019, pp. 1–4.
  * [18] M. Agatonović, Z. Stanković, N. Dončov, L. Sit, B. Milovanović, and T. Zwick, “Application of artificial neural networks for efficient high-resolution 2d doa estimation,” _Radioengineering_ , vol. 21, no. 4, pp. 1178–1186, 2012.
  * [19] Y. Liu, H. Chen, and B. Wang, “Doa estimation based on cnn for underwater acoustic array,” _Applied Acoustics_ , vol. 172, p. 107594, 2021.
  * [20] B. Hu _et al._ , “Doa robust estimation of echo signals based on deep learning networks with multiple type illuminators of opportunity,” _IEEE Access_ , vol. 8, pp. 14 809–14 819, 2020.
  * [21] O. Bialer, N. Garnett, and T. Tirer, “Performance advantages of deep neural networks for angle of arrival estimation,” in _ICASSP 2019-2019 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP)_. IEEE, 2019, pp. 3907–3911.
  * [22] W. Tang, M. Z. Chen, X. Chen, J. Y. Dai, Y. Han, M. Di Renzo, Y. Zeng, S. Jin, Q. Cheng, and T. J. Cui, “Wireless communications with reconfigurable intelligent surface: Path loss modeling and experimental measurement,” _IEEE transactions on wireless communications_ , vol. 20, no. 1, pp. 421–439, 2020.
  * [23] P. Stoica and A. Nehorai, “Music, maximum likelihood, and cramer-rao bound,” _IEEE Transactions on Acoustics, speech, and signal processing_ , vol. 37, no. 5, pp. 720–741, 1989.
  * [24] N. Boumal, _An introduction to optimization on smooth manifolds_. Cambridge University Press, 2023. [Online]. Available: https://www.nicolasboumal.net/book
  * [25] K. Boman and P. Stoica, “Low angle estimation: models, methods, and bounds,” _Digital Signal Processing_ , vol. 11, no. 1, pp. 35–79, 2001.
  * [26] R. Schmidt, “Multiple emitter location and signal parameter estimation,” _IEEE transactions on antennas and propagation_ , vol. 34, no. 3, pp. 276–280, 1986.
  * [27] H. Chen, Y. Wei, Y. Bai, and X. Zhang, “Sparse recovery for doa estimation with a reflection path,” _IEEE Access_ , vol. 6, pp. 70 572–70 581, 2018\.
  * [28] M. Grant and S. Boyd, “Cvx: Matlab software for disciplined convex programming, version 2.1,” 2014.
  * [29] N. Boumal, B. Mishra, P.-A. Absil, and R. Sepulchre, “Manopt, a Matlab toolbox for optimization on manifolds,” _Journal of Machine Learning Research_ , vol. 15, no. 42, pp. 1455–1459, 2014. [Online]. Available: https://www.manopt.org