###### Abstract

Word vector representations are well developed tools for various NLP and
Machine Learning tasks and are known to retain significant semantic and
syntactic structure of languages. But they are prone to carrying and
amplifying bias which can perpetrate discrimination in various applications.
In this work, we explore new simple ways to detect the most stereotypically
gendered words in an embedding and remove the bias from them. We verify how
names are masked carriers of gender bias and then use that as a tool to
attenuate bias in embeddings. Further, we extend this property of names to
show how names can be used to detect other types of bias in the embeddings
such as bias based on race, ethnicity, and age.

Attenuating Bias in Word Vectors

Sunipa Dev Jeff Phillips

University of Utah University of Utah

##  1 BIAS IN WORD VECTORS

Word embeddings are an increasingly popular application of neural networks
wherein enormous text corpora are taken as input and words therein are mapped
to a vector in some high dimensional space. Two commonly used approaches to
implement this are WordToVec [16, 15] and GloVe [17]. These word vector
representations estimate similarity between words based on the context of
their nearby text, or to predict the likelihood of seeing words in the context
of another. Richer properties were discovered such as synonym similarity,
linear word relationships, and analogies such as man : woman :: king : queen.
Their use is now standard in training complex language models.

However, it has been observed that word embeddings are prone to express the
bias inherent in the data it is extracted from [3, 4, 7]. Further, Zhao _et
al._ (2017) [18] and Hendricks _et al._ (2018) [6] show that machine learning
algorithms and their output show more bias than the data they are generated
from.

Word vector embeddings as used in machine learning towards applications which
significantly affect people’s lives, such as to assess credit [11], predict
crime [5], and other emerging domains such judging loan applications and
resumes for jobs or college applications. So it is paramount that efforts are
made to identify and if possible to remove bias inherent in them. Or at least,
we should attempt minimize the propagation of bias within them. For instance,
in using existing word embeddings, Bolukbasi _et al._ (2016) [3] demonstrated
that women and men are associated with different professions, with men
associated with leaderships roles and professions like doctor, programmer and
women closer to professions like receptionist or nurse. Caliskan _et al._
(2017) [7] similarly noted how word embeddings show that women are more
closely associated with arts than math while it is the opposite for men. They
also showed how positive and negative connotations are associated with
European-American versus African-American names.

Our work simplifies, quantifies, and fine-tunes these approaches: we show that
very simple linear projection of all words based on vectors captured by common
names is an effective and general way to significantly reduce bias in word
embeddings. More specifically:

  * 1a.

We demonstrate that simple linear projection of all word vectors along a bias
direction is more effective than the Hard Debiasing of Bolukbasi _et al._
(2016) [3] which is more complex and also partially relies on crowd sourcing.

  * 1b.

We show that these results can be slightly improved by dampening the
projection of words which are far from the projection distance.

  * 2.

We examine the bias inherent in the standard word pairs used for debiasing
based on gender by randomly flipping or swapping these words in the raw text
before creating the embeddings. We show that this alone does not eliminate
bias in word embeddings, corroborating that simple language modification is
not as effective as repairing the word embeddings themselves.

  * 3a.

We show that common names with gender association (e.g., john, amy) often
provides a more effective gender subspace to debias along than using gendered
words (e.g., he, she).

  * 3b.

We demonstrate that names carry other inherent, and sometimes unfavorable,
biases associated with race, nationality, and age, which also corresponds with
bias subspaces in word embeddings. And that it is effective to use common
names to establish these bias directions and remove this bias from word
embeddings.

##  2 DATA AND NOTATIONS

We set as default the text corpus of a Wikipedia dump
(dumps.wikimedia.org/enwiki/latest/enwiki-latest-pages-articles.xml.bz2) with
4.57 billion tokens and we extract a GloVe embedding from it in $D=300$
dimensions per word. We restrict the word vocabulary to the most frequent
$100{,}000$ words. We also modify the text corpus and extract embeddings from
it as described later.

So, for each word in the Vocabulary $W$, we represent the word by the vector
$w_{i}\in\mathbb{R}^{D}$ in the embedding. The bias (e.g., gender) subspace is
denoted by a set of vector $B$. It is typically considered in this work to be
a single unit vector, $v_{B}$ (explained in detail later). As we will revisit,
a single vector is typically sufficient, and will simplify descriptions.
However, these approaches can be generalized to a set of vectors defining a
multi-dimensional subspace.

##  3 HOW TO ATTENUATE BIAS

Given a word embedding, debiasing typically takes as input a set
$\mathcal{E}=\\{E_{1},E_{2},\ldots,E_{m}\\}$ of equality sets. An equality set
$E_{j}$ for instance can be a single pair (e.g., {man, woman}), but could be
more words (e.g., {latina, latino, latinx}) that if the bias connotation (e.g,
gender) is removed, then it would objectively make sense for all of them to be
equal. Our data sets will only use word pairs (as a default the ones in Table
1), and we will describe them as such hereafter for simpler descriptions. In
particular, we will represent each $E_{j}$ as a set of two vectors
$e_{i}^{+},e_{i}^{-}\in\mathbb{R}^{D}$.

Given such a set $\mathcal{E}$ of equality sets, the bias vector $v_{B}$ can
be formed as follows [3]. For each $E_{j}=\\{e_{j}^{+},e_{j}^{-}\\}$ create a
vector $\vec{e}_{i}=e_{i}^{+}-e_{i}^{-}$ between the pairs. Stack these to
form a matrix $Q=[\vec{e}_{1}\;\vec{e}_{2}\;\ldots\;\vec{e}_{m}]$, and let
$v_{B}$ be the top singular vector of $Q$. We revisit how to create such a
bias direction in Section 4.

{man,woman}, {son,daughter}, {he,she}, {his,her},
---
{male,female}, {boy,girl}, {himself,herself},
{guy,gal}, {father,mother}, {john,mary}
Table 1: Gendered Word Pairs

Now given a word vector $w\in W$, we can project it to its component along
this bias direction $v_{B}$ as

$\pi_{B}(w)=\langle w,v_{B}\rangle v_{B}.$

###  3.1 Existing Method : Hard Debiasing

The most notable advance towards debiasing embeddings along the gender
direction has been by Bolukbasi _et al._ (2016) [3] in their algorithm called
Hard Debiasing ($HD$). It takes a set of words desired to be neutralized,
$\\{w_{1},w_{2},\ldots,w_{n}\\}=W_{N}\subset W$, a unit bias subspace vector
$v_{B}$, and a set of equality sets $E_{1},E_{2},\ldots,E_{m}$.

First, words $\\{w_{1},w_{2},\ldots,w_{n}\\}\in W_{N}$ are projected
orthogonal to the bias direction and normalized

$w_{i}^{\prime}=\frac{w_{i}-w_{B}}{||w_{i}-w_{B}||}.$

Second, it corrects the locations of the vectors in the equality sets. Let
$\mu_{j}=\frac{1}{|E|}\sum_{e\in E_{j}}e$ be the mean of an equality set, and
$\mu=\frac{1}{m}\sum_{j=1}^{m}\mu_{j}$ be the mean of of equality set means.
Let $\nu_{j}=\mu-\mu_{j}$ be the offset of a particular equality set from the
mean. Now each $e\in E_{j}$ in each equality set $E_{j}$ is first centered
using their average and then neutralized as

$e^{\prime}=\nu_{j}+\sqrt{1-\|\nu_{j}\|^{2}}\frac{\pi_{B}(e)-v_{B}}{\|\pi_{B}(e)-v_{B}\|}.$

Intuitively $\nu_{j}$ quantifies the amount words in each equality set $E_{j}$
differ from each other in directions apart from the gender direction. This is
used to center the words in each of these sets.

This renders word pairs such as man and woman as equidistant from the neutral
words $w_{i}^{\prime}$ with each word of the pair being centralized and moved
to a position opposite the other in the space. This can filter out properties
either word gained by being used in some other context, like mankind or humans
for the word man.

The word set $W_{N}=\\{w_{1},w_{2},\ldots,w_{n}\\}\subset W$ which is debiased
is obtained in two steps. First it seeds some words as definitionally gendered
via crowd sourcing and using dictionary definitions; the complement – ones not
selected in this step – are set as neutral. Next, using this seeding an SVM is
trained and used to predict among all $W$ the set of other biased $W_{B}$ or
neutral words $W_{N}$. This set $W_{N}$ is taken as desired to be neutral and
is debiased. Thus not all words $W$ in the vocabulary are debiased in this
procedure, only a select set chosen via crowd-sourcing and definitions, and
its extrapolation. Also the word vectors in the equality sets are also handled
separately. This makes this approach not a fully automatic way to debias the
vector embedding.

###  3.2 Alternate and Simple Methods

We next present some simple alternatives to HD which are simple and fully
automatic. These all assume a bias direction $v_{B}$.

Subtraction.  As a simple baseline, for _all_ word vectors $w$ subtract the
gender direction $v_{B}$ from $w$:

$w^{\prime}=w-v_{B}.$

Linear Projection.  A better baseline is to project _all_ words $w\in W$
orthogonally to the bias vector $v_{B}$.

$w^{\prime}=w-\pi_{B}(w)=w-\langle w,v_{B}\rangle v_{B}.$

This enforces that the updated set $W^{\prime}=\\{w^{\prime}\mid w\in W\\}$
has no component along $v_{B}$, and hence the resulting span is only $D-1$
dimensions. Reducing the total dimension from say $300$ to $299$ should have
minimal effects of expressiveness or generalizability of the word vector
embeddings.

Bolukbasi _et al._ [3] apply this same step to a dictionary definition based
extrapolation and crowd-source-chosen set of word pairs $W_{N}\subset W$. We
quantify in Section 5 that this single universal projection step debiases
better than HD.

For example, consider the bias as gender, and the equality set with words man
and woman. Linear projection will subtract from their word embeddings the
proportion that were along the gender direction $v_{B}$ learned from a larger
set of equality pairs. It will make them close-by but not exactly equal. The
word man is used in many extra senses than the word woman; it is used to refer
to humankind, to a person in general, and in expressions like “oh man”. In
contrast a simpler word pair with fewer word senses, like (he \- she) and (him
\- her), we can expect them to be almost at identical positions in the vector
space after debiasing, implying their synonymity.

Thus, this approach uniformly reduces the component of the word along the bias
direction without compromising on the differences that words (and word pairs)
have.

###  3.3 Partial Projection

A potential issue with the simple approaches is that they can significantly
change some embedded words which are definitionally biased (e.g., the neutral
words $W_{B}$ described by Bolukbasi _et al._ [3]). _[[We note that this may
not *actually* be a problem (see Section 5); the change may only be associated
with the bias, so removing it would then not change the meaning of those words
in any way except the ones we want to avoid.]]_ However, these intuitively
should be words which have correlation with the bias vector, but also are far
in the orthogonal direction. In this section we explore how to automatically
attenuate the effect of the projection on these words.

This stems from the observation that given a bias direction, the words which
are most extreme in this direction (have the largest dot product) sometimes
have a reasonable biased context, but some do not. These “false positives” may
be large normed vectors which also happen to have a component in the bias
direction.

Figure 1: Illustration of $\eta$ and $\beta$ for word vector $w$.

We start with a bias direction $v_{B}$ and mean $\mu$ derived from equality
pairs (defined the same way as in context of HD). Now given a word vector $w$
we decompose it into key values along two components, illustrated in Figure 1.
First, we write its bias component as

$\beta(w)=\langle w,v_{B}\rangle-\langle\mu,v_{B}\rangle.$

This is the difference of $w$ from $\mu$ when both are projected onto the bias
direction $v_{B}$.

Second, we write a (residual) orthogonal component

$r(w)=w-\langle w,v_{B}\rangle v_{B}.$

Let $\eta(w)=\|r(w)\|$ be its value. It is the orthogonal distance from the
bias vector $v_{B}$; recall we chose $v_{B}$ to pass through the origin, so
the choice of $\mu$ does not affect this distance.

Now we will maintain the orthogonal component ($r(w)$, which is in a subspace
spanned by $D-1$ out of $D$ dimensions) but adjust the bias component
$\beta(w)$ to make it closer to $\mu$. But the adjustment will depend on the
magnitude $\eta(w)$. As a default we set

$w^{\prime}=\mu+r(w)$

so all word vectors retain their orthogonal component, but have a fixed and
constant bias term. This is functionally equivalent to the Linear Projection
approach; the only difference is that instead of having a $0$ magnitude along
$v_{B}$ (and the orthogonal part unchanged), it instead has a magnitude of
constant $\mu$ along $v_{B}$ (and the orthogonal part still unchanged). This
adds a constant to every inner product, and a constant offset to any linear
projection or classifier. If we are required to work with normalize vectors
(we do not recommend this as the vector length captures veracity information
about its embedding), we can simple set $w^{\prime}=r(w)/\|r(w)\|$.

Given this set-up, we now propose three modifications. In each set

$w^{\prime}=\mu+r(w)+\beta\cdot f_{i}(\eta(w))\cdot v_{B}$

were $f_{i}$ for $i=\\{1,2,3\\}$ is a function of only the orthogonal value
$\eta(w)$. For the default case $f(\eta)=0$

$\displaystyle f_{1}(\eta)$ $\displaystyle=\sigma^{2}/(\eta+1)^{2}$
$\displaystyle f_{2}(\eta)$ $\displaystyle=\exp(-\eta^{2}/\sigma^{2})$
$\displaystyle f_{3}(\eta)$ $\displaystyle=\max(0,\sigma/2\eta)$

Here $\sigma$ is a hyperparameter that controls the importance of $\eta$; in
Section 3.4 we show that we can just set $\sigma=1$. x

Figure 2: The gendered region as per the three variations of projection. Both
points P1 and P2 have a dot product of 1.0 initially with the gender subspace.
But their orthogonal distance to it differs, as expressed by their dot product
with the other 299 dimensions.

In Figure 2 we see the regions of the $(\eta,\beta)$-space that the functions
$f$, $f_{1}$ and $f_{2}$ consider gendered. $f$ projects all points onto the
$y=\mu$ line. But variants $f_{1}$, $f_{2}$, and $f_{3}$ are represented by
curves that dampen the bias reduction to different degrees as $\eta$
increases. Points P1 and P2 have the same dot products with the bias direction
but different dot products along the other $D-1$ dimensions. We can observe
the effects of each dampening function as $\eta$ increases from P1 to P2.

###  3.4 SETTING $\sigma=1$.

To complete the damping functions $f_{1}$, $f_{2}$, and $f_{3}$, we need a
value $\sigma$. If $\sigma$ is larger, then more word vectors have bias
completely removed; if $\sigma$ is smaller, than more words are unaffected by
the projection. The goal is that words $S$ which are likely to carry a bias
connotation to have little damping (small $f_{i}$ values) and words $T$ which
are unlikely to carry a bias connotation to have more damping (large $f_{i}$
values – they are not moved much).

Given sets $S$ and $T$, we can define a gain function

$\gamma_{i,\rho}(\sigma)=\sum_{s\in
S}\beta(s)(1-f_{i}(\eta(s)))-\rho\sum_{t\in T}\beta(t)(1-f_{i}(\eta(t))),$

with a regularization term $\rho$. The gain $\gamma$ is large when most bias
words in $S$ have very little damping (small $f_{i}$, large $1-f_{i}$), and
the opposite is true for the neutral words in $T$. We want the neutral words
to have large $f_{i}$ and hence small $1-f_{i}$, so they do not change much.

To define the gain function, we need sets $S$ and $T$; we do so with the bias
of interest as gender. The biased set $S$ is chosen among a set of $1000$
popular names in $W$ which (based on babynamewizard.com and and SSN databases
[2, 1]) are strongly associated with a gender. The neutral set $T$ is chosen
as the most frequent $1000$ words from $W$, after filtering out obviously
gendered words like names man and he. We also omit occupation words like
doctor and others which may carry unintentional gender bias (these are words
we would like to automatically de-bias). The neutral set may not be perfectly
non-gendered, but it provides reasonable approximation of all non-gendered
words.

We find for an array of choices for $\rho$ (we tried $\rho=1$, $\rho=10$, and
$\rho=100$), the value $\sigma=1$ approximately maximizes the gain function
$\gamma_{i,\rho}(\sigma)$ for each $i\in\\{1,2,3\\}$. So for hereafter we fix
$\sigma=1$.

Although these sets $S$ and $T$ play a role somewhat similar to the crowd-
sourced sets $W_{B}$ and $W_{N}$ from HD that we hoped to avoid, the role here
is much reduced. This is just to verify that a choice of $\sigma=1$ is
reasonable, and otherwise they are not used.

123 $0.00$ $0.50$ $0.75$ $1.00$

Figure 3:  Fractional singular values for avg male - female words (as per
Table 1) after flipping with probability (from left to right) $0.0$ (the
original data set), $0.5$, $0.75$, and $1.0$.

###  3.5 Flipping the Raw Text

Since the embeddings preserve inner products of the data from which it is
drawn, we explore if we can make the data itself gender unbiased and then
observe how that change shows up in the embedding. Unbiasing a textual corpus
completely can be very intricate and complicated since there are a many
(sometimes implicit) gender indicators in text. Nonetheless, we propose a
simple way of neutralizing bias in textual data by using word pairs
$E_{1},E_{2},\ldots E_{m}$; in particular, when we observe in raw text on part
of a word part, we randomly flip it to the other pair. For instance for
gendered word pairs (e.g., (he \- she)) in a string “he was a doctor” we may
flip to “she was a doctor.”

We implement this procedure over the entire input raw text, and try various
probabilities of flipping each observed word, focusing on probabilities $0.5$,
$0.75$ and $1.00$. The first $0.5$-flip probability makes each element of a
word pair equally likely. The last $1.00$-flip probability reverses the roles
of those word pairs, and $0.75$-flip probability does something in between. We
perform this set of experiments on the default Wikipedia data set and switch
between word pairs (say man $\rightarrow$ woman, she $\rightarrow$ he, etc),
from a list larger that Table 3 consisting of 75 word pairs; see Supplementary
Material D.1.

We observe how the proportion along the principal component changes with this
flipping in Figure 3. We see that flipping with $0.5$ somewhat dampens the
difference between the different principal components. On the other hand
flipping with probability $1.0$ (and to a lesser extent $0.75$) exacerbates
the gender components rather than dampening it. Now there are two components
significantly larger than the others. This indicates this flipping is only
addressing part of the explicit bias, but missing some implicit bias, and
these effects are now muddled.

We list some gender biased analogies in the default embedding and how they
change with each of the methods described in this section in Table 2.

Table 2: What analogies look like before and after damping gender by different methods discussed : hard Debiaisng, flipping words in text corpus, subtraction and projection analogy head | original | HD | flipping | subtraction | projection
---|---|---|---|---|---
|  |  | 0.5 | 0.75 | 1.0 |  |
man : woman :: doctor : | nurse | surgeon | dr | dr | medicine | physician | physician
man : woman :: footballer : | politician | striker | midfielder | goalkeeper | striker | politician | midfielder
he : she :: strong : | weak | stronger | weak | strongly | many | well | stronger
he : she :: captain : | mrs | lieutenant | lieutenant | colonel | colonel | lieutenant | lieutenant
john : mary :: doctor : | nurse | physician | medicine | surgeon | nurse | father | physician

##  4 THE BIAS SUBSPACE

We explore ways of detecting and defining the bias subspace $v_{B}$ and
recovering the most gendered words in the embedding. Recall as default, we use
$v_{B}$ as the top singular vector of the matrix defined by stacking vectors
$\vec{e}_{i}=e_{i}^{+}-e_{i}^{-}$ of biased word pairs. We primarily focus on
gendered bias, using words in Table 1, and show later how to effectively
extend to other biases. We discuss this in detail in Supplementary Material C.

#### Most gendered words.

The dot product, $\langle v_{B},w\rangle$ of the word vectors $w$ with the
gender subspace $v_{B}$ is a good indicator of how gendered a word is. The
magnitude of the dot product tells us of the length along the gender subspace
and the sign tells us whether it is more female or male. Some of the words
denoted as most gendered are listed in Table 3.

Table 3: Some of the most gendered words in default embedding; and most
gendered adjectives and occupation words.

Gendered Words
---
miss | herself | forefather | himself
maid | heroine | nephew | congressman
motherhood | jessica | zahir | suceeded
adriana | seductive | him | sir

Female Adjectives | Male Adjectives
---|---
glamorous | strong
diva | muscular
shimmery | powerful
beautiful | fast

Female Occupations | Male Occupations
---|---
nurse | soldier
maid | captain
housewife | officer
prostitute | footballer

###  4.1 Bias Direction using Names

When listing gendered words by $|\langle v_{B},w\rangle|$, we observe that
many gendered words are names. This indicates the potential to use names as an
alternative (and potentially in a more general way) to bootstrap finding the
gender direction.

From the top 100K words, we extract the 10 most common male
$\\{m_{1},m_{2},\ldots,m_{10}\\}$ and female $\\{s_{1},s_{2},\ldots,s_{10}\\}$
names which are not used in ambiguous ways (e.g., not the name hope which
could also refer to the sentiment). We pair these 10 names from each category
(male, female) randomly and compute the SVD as before. We observe in 4 that
the fractional singular values show a similar pattern as with the list of
correctly gendered word pairs like (man \- woman), (he \- she), etc.

Figure 4: Proportion of singular values along principal directions (left)
using names as indicators, and (right) using word pairs from Table 1 as
indicators

But this way of pairing names is quite imprecise. These names are not
‘opposites’ of each other in the sense that word pairs are. So, we modify how
we compute $v_{B}$ now so that we can better use names to detect the bias in
the embedding. The following method gives us this advantage where we do not
necessarily need word pairs or equality sets as in Bolukbasi _et al._ [3].

Our gender direction is calculated as,

$v_{B,\textsf{names}}=\frac{s-m}{\|s-m\|},$

where $s=\frac{1}{10}\sum_{i}s_{i}$ and $m=\frac{1}{10}\sum_{i}m_{i}$.

Using the default Wikipedia dataset, we found that this is a good approximator
of the gender subspace defined by the first right singular vector calculated
using gendered words from Table 1; there dot product is $0.809$. We find
similar large dot product scores for other datasets too.

Here too we collect all the most gendered words as per the gender direction
$v_{B,\textsf{names}}$ determined by these names. Most gendered words returned
are similar as using the default $v_{B}$, like occupational words, adjectives,
and synonyms for each gender. We find names to express similar classification
of words along male - female vectors with homemaker more female and policeman
being more male. We illustrate this in more detail in Table 4.

Table 4: Gendered occupations as observed in word embeddings using names as the gender direction indicator  Female Occ | Male Occ | Female* Occ | Male* Occ
---|---|---|---
nurse | captain | policeman | policeman
maid | cop | detective | cop
actress | boss | character | character
housewife | officer | cop | assassin
dancer | actor | assassin | bodyguard
nun | scientist | actor | waiter
waitress | gangster | waiter | actor
scientist | trucker | butler | detective

Using that direction, we debias by linear projection. There is a similar shift
in analogy results. We see a few examples in Table 5.

Table 5: What analogies look like before and after removing the gender direction using names analogy head | original | subtraction | projection
---|---|---|---
man : woman :: doctor : | nurse | physician | physician
man : woman :: footballer : | politician | politician | midfielder
he : she :: strong : | weak | very | stronger
he : she :: captain : | mrs | lieutenant | lieutenant
john : mary :: doctor : | nurse | dr | dr

##  5 QUANTIFYING BIAS

In this section we develop new measures to quantify how much bias has been
removed from an embedding, and evaluate the various techniques we have
developed for doing so.

As one measure, we use the Word Embedding Association Test (WEAT) test
developed by Caliskan _et al._ (2017) [7] as analogous to the IAT tests to
evaluate the association of male and female gendered words with two categories
of target words: career oriented words versus family oriented words. We detail
WEAT and list the exact words used (as in [7]) in Supplementary Material B;
smaller values are better.

Bolukbasi _et al._ [3] evaluated embedding bias use a crowsourced judgement of
whether an analogy produced by an embedding is biased or not. Our goal was to
avoid crowd sourcing, so we propose two more automatic tests to qualitatively
and uniformly evaluate an embedding for the presence of gender bias.

#### Embedding Coherence Test (ECT).

A way to evaluate how the neutralization technique affects the embedding is to
evaluate how the nearest neighbors change for (a) gendered pairs of words
$\mathcal{E}$ and (b) indirect-bias-affected words such as those associated
with sports or occupational words (e.g., football, captain, doctor). We use
the gendered word pairs in Table 1 for $\mathcal{E}$ and the professions list
$P=\\{p_{1},p_{2},\ldots,p_{k}\\}$ as proposed and used by Bolukbasi _et al._
https://github.com/tolga-b/debiaswe (see also Supplementary Material D.2) to
represent (b).

  * S1:

For all word pair $\\{e_{j}^{+},e_{j}^{-}\\}=E_{j}\in\mathcal{E}$ we compute
two means $m=\frac{1}{|\mathcal{E}|}\sum_{E_{j}\in\mathcal{E}}e_{j}^{+}$ and
$s=\frac{1}{|\mathcal{E}|}\sum_{E_{j}\in\mathcal{E}}e_{j}^{-}$. We find the
cosine similarity of both $m$ and $s$ to all words $p_{i}\in P$. This creates
two vectors $u_{m},u_{s}\in\mathbb{R}^{k}$.

  * S2:

We transform these similarity vectors to replace each coordinate by its rank
order, and compute the Spearman Coefficient (in $[-1,1]$, larger is better)
between the rank order of the similarities to words in $P$.

Thus, here, we care about the order in which the words in $P$ occur as
neighbors to each word pair rather than the exact distance. The exact distance
between each word pair would depend on the usage of each word and thus on all
the different dimensions other than the gender subspace too. But the order
being relatively the same, as determined using Spearman Coefficient would
indicate the dampening of bias in the gender direction (i.e., if doctor by
profession is the 2nd closest of all professions to both man and woman, then
the embedding has a dampened bias for the word doctor in the gender
direction). Neutralization should ideally bring the Spearman coefficient
towards 1.

#### Embedding Quality Test (EQT).

The demonstration by Bolukbasi _et al._ [3] about the skewed gender roles in
embeddings using analogies is what we try to quantify in this test. We attempt
to quantify the improvement in analogies with respect to bias in the
embeddings.

We use the same sets $\mathcal{E}$ and $P$ as in the ECT test. However, for
each profession $p_{i}\in P$ we create a list $S_{i}$ of their plurals and
synonyms from WordNet on NLTK [14].

  * S1:

For each word pair $\\{e_{j}^{+},e_{j}^{-}\\}=E_{j}\in\mathcal{E}$, and each
occupation word $p_{i}\in P$, we test if the analogy
$e_{j}^{+}:e_{j}^{-}::p_{i}$ returns a word from $S_{i}$. If yes, we set
$Q(E_{j},p_{i})=1$, and $Q(E_{j},p_{i})=0$ otherwise.

  * S2:

Return the average value across all combinations
$\frac{1}{|\mathcal{E}|}\frac{1}{k}\sum_{E_{j}\in\mathcal{E}}\sum_{p_{i}\in
P}Q(E_{j},p_{i})$.

The scores for EQT are typically much smaller than for ECT. We explain two
reasons for this.

First, EQT does not check for if the analogy makes relative sense, biased or
otherwise. So, “man : woman :: doctor : nurse” is as wrong as “man : woman ::
doctor : chair.” This pushes the score down.

Second, synonyms in each set $s_{i}$ as returned by WordNet [8] on the Natural
Language Toolkit, NLTK [14] do not always contain all possible variants of the
word. For example, the words psychiatrist and psychologist can be seen as
analogous for our purposes here but linguistically are removed enough that
WordNet does not put them as synonyms together. Hence, even after debiasing,
if the analogy returns “man : woman :: psychiatrist : psychologist“ S1 returns
0. Further, since the data also has several misspelt words, archeologist is
not recognized as a synonym or alternative for the word archaeologist. For
this too S1 returns a 0.

The first caveat can be side-stepped by restricting the pool of words we
search over for the analogous word to be from list $P$. But it is debatable if
an embedding should be penalized equally for returning both nurse or chair for
the analogy “man : woman :: doctor : ?”

This measures the quality of analogies, with better quality having a score
closer to $1$.

Table 6: Performance on ECT, EQT and WEAT by the different debiasing methods; and performance on standard similarity and analogy tests. analogy head | original | HD | flipping | subtraction | projection
---|---|---|---|---|---
|  |  | 0.5 | 0.75 | 1.0 | word pairs | names | word pairs | names
ECT (word pairs) | 0.798 | 0.917 | 0.983 | 0.984 | 0.683 | 0.963 | 0.936 | 0.996 | 0.943
ECT (names) | 0.832 | 0.968 | 0.714 | 0.662 | 0.587 | 0.923 | 0.966 | 0.935 | 0.999
EQT | 0.128 | 0.145 | 0.131 | 0.098 | 0.085 | 0.268 | 0.236 | 0.283 | 0.291
WEAT | 1.623 | 1.221 | 1.164 | 1.09 | 1.03 | 1.427 | 1.440 | 1.233 | 1.219
WSim | 0.637 | 0.537 | 0.567 | 0.537 | 0.536 | 0.627 | 0.636 | 0.627 | 0.629
Simlex | 0.324 | 0.314 | 0.317 | 0.314 | 0.264 | 0.302 | 0.312 | 0.321 | 0.321
Google Analogy | 0.623 | 0.561 | 0.565 | 0.561 | 0.321 | 0.538 | 0.565 | 0.565 | 0.584

#### Evaluating embeddings.

We mainly run $4$ methods to evaluate our methods WEAT, EQT, and two variants
of ECT: ECT (word pairs) uses $\mathcal{E}$ defined by words in Table 1 and
ECT (names) which uses vectors $m$ and $s$ derived by gendered names.

We observe in Table 6 that the ECT score increases for all methods in
comparison to the non-debiased (the original) word embedding; the exception is
flipping with $1.0$ probability score for ECT (word pairs) and all flipping
variants for ECT (names). Flipping does nothing to affect the names, so it is
not surprising that it does not improve this score; further indicating that it
is challenging to directly fix bias in raw text before creating embeddings.
Moreover, HD has the lowest score (of $0.917$) whereas projection obtains
scores of $0.996$ (with $v_{B}$) and $0.943$ (with $v_{B,\textsf{names}}$).

EQT is a more challenging test, and the original embedding only achieves a
score of $0.128$, and HD only obtains $0.145$ (that is $12-15\%$ of occupation
words have their related word as nearest neighbor). On the other hand,
projection increases this percentage to $28.3\%$ (using $v_{B}$) and $29.1\%$
(using $v_{B,\textsf{names}}$). Even subtraction does nearly as well at
between $23-27\%$. Generally, the subtraction always performs slightly worse
than projection.

For the WEAT test, the original data has a score of $1.623$, and this is
decreased the most by all forms of flipping, down to about $1.1$. HD and
projection do about the same with HD obtaining a score of $1.221$ and
projection obtaining $1.219$ (with $v_{B,\textsf{names}}$) and $1.234$ (with
$v_{B}$); values closer to 0 are better (See Supplementary Material B).

In the bottom of Table 6 we also run these approaches on standard similarity
and analogy tests for evaluating the quality of embeddings. We use cosine
similarity [13] on WordSimilarity-353 (WSim, 353 word pairs) [9] and
SimLex-999 (Simlex, 999 word pairs) [10], each of which evaluates a Spearman
coefficient (larger is better). We also use the Google Analogy Dataset using
the function 3COSADD [12] which takes in three words which for a part of the
analogy and returns the 4th word which fits the analogy the best.

We observe (as expected) that all that debiasing approaches reduce these
scores. The largest decrease in scores (between $1\%$ and $10\%$) is almost
always from HD. Flipping at $0.5$ rate is comparable to HD. And simple linear
projection decreases the least (usually only about $1\%$, except on analogies
where it is $7\%$ (with $v_{B}$) or $5\%$ (with $v_{B,\textsf{names}}$).

In Table 7 we also evaluate the damping mechanisms defined by $f_{1}$,
$f_{2}$, and $f_{3}$, using $v_{B}$. These are very comparable to simple
linear projection (represented by $f$). The scores for ECT, EQT, and WEAT are
all about the same as simple linear projection, usually slightly worse.

Table 7:  Performance of damped linear projection using word pairs.  Tests | f | $f_{1}$ | $f_{2}$ | $f_{3}$
---|---|---|---|---
ECT | 0.996 | 0.994 | 0.995 | 0.997
EQT | 0.283 | 0.280 | 0.292 | 0.287
WEAT | 1.233 | 1.253 | 1.245 | 1.241
WSim | 0.627 | 0.628 | 0.627 | 0.627
Simlex | 0.321 | 0.324 | 0.324 | 0.324
Google Analogy | 0.565 | 0.571 | 0.569 | 0.569

While ECT, EQT and WEAT scores are in a similar range for all of $f$, $f_{1}$,
$f_{2}$, and $f_{3}$; the dampened approaches $f_{1}$, $f_{2}$, and $f_{3}$
performs better on the Google Analogy test. This test set is devoid of bias
and is made up of syntactic and semantic analogies. So, a score closer to that
of the original, biased embedding, tells us that more structure has been
retained by $f_{1}$, $f_{2}$ and $f_{3}$. Overall, any of these approaches
could be used if a user wants to debias while retaining as much structure as
possible, but otherwise linear projection (or $f$) is roughly as good as these
dampened approaches.

##  6 DETECTING OTHER BIAS USING NAMES

We saw so far how projection combined with finding the gender direction using
names works well and works as well as projection combined with finding the
gender direction using word pairs.

We explore here a way of extending this approach to detect other kinds of bias
where we cannot necessarily find good word pairs to indicate a direction, like
Table 1 for gender, but where names are known to belong to certain protected
demographic groups. For example, there is a divide between names that
different racial groups tend to use more. Caliskan _et al._ [7] use a list of
names that are more African-American (AA) versus names that are more European-
American (EA) for their analysis of bias. There are similar lists of names
that are distinctly and commonly used by different ethnic, racial (e.g.,
Asian, African-American) and even religious (for e.g., Islamic) groups.

We first try this with two common demographic group divides : Hispanic /
European-American and African-American / European-American.

Hispanic and European-American names. Even though we begin with most commonly
used Hispanic (H) names (Supplementary Material D), this is tricky as not all
names occur as much as European American names and are thus not as well
embedded. We use the frequencies from the dataset to guide us in selecting
commonly used names that are also most frequent in the Wikipedia dataset.
Using the same method as Section 4.1, we determine the direction,
$v_{B,\textsf{names}}$, which encodes this racial difference and find the
words most commonly aligned with it. Other Hispanic and European-American
names are the closest words. But other words like, latino or hispanic also
appear to be close, which affirms that we are capturing the right subspace.

African-American and European-American names. We see a similar trend when we
use African-American names and European-American names (Figure 5). We use the
African-American names used by Caliskan _et al._ (2017) [7]. We determine the
bias direction by using method in Section 4.1.

We plot in Figure 5 a few occupation words along the axes defined by H-EA and
AA-EA bias directions, and compare them with those along the male-female axis.
The embedding is different among the groups, and likely still generally more
subordinate-biased towards Hispanic and African-American names as it was for
female. Although footballer is more Hispanic than European-American, while
maid is more neutral in the racial bias setting than the gender setting. We
see this pattern repeated across embeddings and datasets (see Supplementary
Material A).

When we switch the type of bias, we also end up finding different patterns in
the embeddings. In the case of both of these racial directions, there is a the
split in not just occupation words but other words that are detected as highly
associated with the bias subspace. It shows up foremost among the closest
words of the subspace of the bias. Here, we find words like drugs and illegal
close to the H-EA direction while, close to the AA-EA direction, we retrieve
several slang words used to refer to African-Americans. These word
associations with each racial group can be detected by the WEAT tests (lower
means less bias) using positive and negative words as demonstrated by Caliskan
_et al._ (2017) [7]. We evaluate using the WEAT test before and after linear
projection debiasing in Table 8. For each of these tests, we use half of the
names in each category for finding the bias direction and the other half for
WEAT testing. This selection is done arbitrarily and the scores are averaged
over 3 such selections.

Figure 5:  Gender and racial bias in the embedding Table 8:  WEAT positive-negative test scores before and after debiasing | Before Debiasing | After Debiasing
---|---|---
EA-AA | 1.803 | 0.425
EA-H | 1.461 | 0.480
Youth-Aged | 0.915 | 0.704

More qualitatively, as a result of the dampening of bias, we see that biased
words like other names belonging to these specific demographic groups, slang
words, colloquial terms like latinos are removed from the closest 10% words.
This is beneficial since the distinguishability of demographic characteristics
based on names is what shows up in these different ways like in occupational
or financial bias.

#### Age-associated names.

We observed that names can be masked carriers of age too. Using the database
for names through time [1] and extracting the most common names from early
1900s as compared to late 1900s and early 2000s, we find a correlation between
these names (see Supplementary Material) and age related words. In Figure 6,
we see a clear correlation between age and names. Bias in this case does not
show up in professions as clearly as in gender but in terms of association
with positive and negative words [7]. We again evaluate using a WEAT test in
Table 8, the bias before and after debiasing the embedding.

Figure 6: Detecting Age with Names : a plot of age related terms along names
from different centuries.

##  7 DISCUSSION

Different types of bias exist in textual data. Some are easier to detect and
evaluate. Some are harder to find suitable and frequent indicators for and
thus, to dampen. Gendered word pairs and gendered names are frequent enough in
textual data to allow us to successfully measure it in different ways and
project the word embeddings away from the subspace occupied by gender. Other
types of bias don’t always have a list of word pairs to fall back on to do the
same. But using names, as we see here, we can measure and detect the different
biases anyway and then project the embedding away from. In this work we also
see how a weighted variant of propection removes bias while retaining best the
inherent structure of the word embedding.

## References

  * [1] https://www.ssa.gov/oact/babynames/.
  * [2] http://www.babynamewizard.com.
  * [3] T Bolukbasi, K W Chang, J Zou, V Saligrama, and A Kalai.  Man is to computer programmer as woman is to homemaker? debiasing word embeddings.  In ACM Transactions of Information Systems, 2016.
  * [4] T Bolukbasi, K W Chang, J Zou, V Saligrama, and A Kalai.  Quantifying and reducing bias in word embeddings.  2016\.
  * [5] Tim Brennan, William Dieterich, and Beate Ehret.  Evaluating the predictive validity of the compas risk and needs assessment system.  Criminal Justice and Behavior, 36(1):21–40, 2009.
  * [6] Kaylee Burns, Lisa Anne Hendricks, Trevor Darrell, and Anna Rohrbach.  Women also snowboard: Overcoming bias in captioning models.  CoRR, abs/1803.09797, 2018.
  * [7] Aylin Caliskan, Joanna J. Bryson, and Arvind Narayanan.  Semantics derived automatically from language corpora contain human-like biases.  Science, 356(6334):183–186, 2017.
  * [8] Christiane Fellbaum.  WordNet: An Electronic Lexical Database.  Bradford Books, 1998.
  * [9] L Finkelstein, E Gabrilovich, Y Matias, E Rivlin, Z Solan, G Wolfman, and etal.  Placing search in context : The concept revisited.  In ACM Transactions of Information Systems, volume 20, pages 116–131, 2002.
  * [10] F Hill, R Reichart, and A Korhonen.  Simlex-999 : Evaluating semantic models with (genuine) similarity estimation.  In Computational Linguistics, volume 41, pages 665–695, 2015.
  * [11] Amir E. Khandani, Adlar J. Kim, and Andrew Lo.  Consumer credit-risk models via machine-learning algorithms.  Journal of Banking & Finance, 34(11):2767–2787, 2010.
  * [12] Omer Levy and Yoav Goldberg.  Neural word embedding as implicit matrix factorization.  In NIPS, 2013.
  * [13] Omer Levy and Yoav Goldberg.  Linguistic regularities of sparse and explicit word representations.  In CoNLL, 2014.
  * [14] Edward Loper and Steven Bird.  Nltk: The natural language toolkit.  In Proceedings of the ACL-02 Workshop on Effective Tools and Methodologies for Teaching Natural Language Processing and Computational Linguistics - Volume 1, ETMTNLP ’02, pages 63–70, Stroudsburg, PA, USA, 2002\. Association for Computational Linguistics.
  * [15] Tomas Mikolov, Kai Chen, Greg Corrado, and Jeffrey Dean.  Efficient estimation of word representations in vector space.  Technical report, arXiv:1301.3781, 2013.
  * [16] Tomas Mikolov, Ilya Sutskever, Kai Chen, Greg Corrado, and Jeffrey Dean.  Distributed representations of words and phrases and their compositionality.  In NIPS, pages 3111–3119, 2013.
  * [17] Jeffrey Pennington, Richard Socher, and Christopher D. Manning.  Glove: Global vectors for word representation.  2014\.
  * [18] Jieyu Zhao, Tianlu Wang, Mark Yatskar, Vicente Ordonez, and Kai-Wei Chang.  Men also like shopping: Reducing gender bias amplification using corpus-level constraints.  CoRR, abs/1707.09457, 2017.

Supplementary material for:
Attenuating Bias in Word Embeddings

##  Appendix A Bias in different embeddings

We explore here how gender bias is expressed across different embeddings,
datasets and embedding mechanisms. Similar patterns are reflected across all
as seen in Figure 7.

For this verification of the permeative nature of bias across datasets and
embeddings, we use the GloVe embeddings of a Wikipedia dump
(dumps.wikimedia.org/enwiki/latest/enwiki-latest-pages-articles.xml.bz2,
$4.7$B tokens) Common Crawl ($840$B tokens, $2.2$M vocab) and Twitter ($27$B
tokens, $1.2$M vocab) from https://nlp.stanford.edu/projects/glove/, and the
WordToVec embedding of Google News ($100$B tokens, $3$M vocab) from
https://code.google.com/archive/p/word2vec/.

(a)
(b)
(c)

(d)

Figure 7: Gender Bias in Different Embeddings : (a) GloVeon Wikipedia, (b)
GloVeon Twitter, (c) GloVeon Common Crawl and (d) WordToVecon GoogleNews
datasets.

(a)  (b)  (c)  (d)

Figure 8: Bias in adjectives along the gender direction : (a) GloVe on default
Wikipedia dataset, (b) GloVe on Common Crawl (840B token dataset), (c) GloVe
on Twitter dataset and (d) WordToVec on Google News dataset

(a)  (b)
(c)  (d)

Figure 9: Racial bias in different embeddings : Occupation words along the
European American - Hispanic axis in GloVe embeddings of (a) Common Crawl and
(b) Twitter Dataset and the European American - African American axis in GloVe
embeddings of (a) Common Crawl and (b) Twitter Dataset

##  Appendix B WORD EMBEDDING ASSOCIATION TEST

Word Embedding Association Test (WEAT) was defined as an analogue to Implicit
Association Test (IAT) by Caliskan _et al._ [7]. It checks for human like bias
associated with words in word embeddings. For example, it found career
oriented words (executive, career, etc) more associated with male names and
male gendered words (’man’,’boy’ etc) than female names and gendered words and
family oriented words (’family’,’home’ etc) more associated with female names
and words than male. We list a set of words used for WEAT by Calisan _et al._
and that we used in our work below.

For two sets of target words X and Y and attribute words A and B, the WEAT
test statistic is :

$s(X,Y,A,B)=\sum_{x\in X}s(x,A,B)-s(y,A,B)$

where,

$s(w,A,B)=mean_{a\in A}cos(a,w)-mean_{b\in B}cos(b,w)$ and, $cos(a,b)$ is the
cosine distance between vector a and b.

This score is normalized by $std-dev_{w\in X\cup Y}s(w,A,B)$. So, closer to 0
this value is, the less bias or preferential association target word groups
have to the attribute word groups.

Here target words are occupation words or career/family oriented words and
attributes are male/female words or names.

Career : { executive, management, professional, corporation, salary, office,
business, career }

Family : { home, parents, children, family, cousins, marriage, wedding,
relatives }

Male names : { john, paul, mike, kevin, steve, greg, jeff, bill }

Female names : { amy, joan, lisa, sarah, diana, kate, ann, donna }

Male words : { male, man, boy, brother, he, him, his, son }

Female words : { female, woman, girl, she, her, hers, daughter }

##  Appendix C DETECTING THE GENDER DIRECTION

For this, we take a set of gendered word pairs as listed in Table 3. From our
default Wikipedia dataset, using the embedded vectors for these word pairs
(i.e., (woman \- man), (she \- he), etc), we create a basis for the subspace
$F$, of dimension $10$. We then try to understand the distribution of variance
in this subspace. To do so, we project the entire dataset onto this subspace
$F$, and take the SVD. The top chart in Figure 10 shows the singular values of
the entire data in this subspace $F$. We observe that there is a dominant
first singular vector/value which is almost twice the size of the second
value. After the this drop, the decay is significantly more gradual. This
suggests to use only the top singular vector of $F$ as the gender subspace,
not $2$ or more of these vectors.

To grasp how much of this variation is from the correlation along the gender
direction, and how much is just random variation, we repeat this experiment,
again in Figure 10, with different ways of creating the subspace $F$. First in
chart (b), we generate $10$ vectors, with one word chosen randomly, and one
chosen from the gendered set (e.g., chair-woman). Second in chart (c), we
generate $10$ vectors between two random words from our set of the $100{,}000$
most frequent words; these are averaged over $100$ random iterations due to
higher variance in the plots. Finally in chart (d), we generate $10$ random
unit vectors in $\mathbb{R}^{300}$. We observe that the pairs with one
gendered vector in each pair still exhibits a significant drop in singular
values, but not as drastic as with both pairs. The other two approaches have
no significant drop since they do not in general contain a gendered word with
interesting subspace. All remaining singular values, and their decay appears
similar to the non-leading ones from the charts (a) and (b). This further
indicates that there is roughly one important gender direction, and any
related subspace is not significantly different than a random one in the word-
vector embedding.

Now, for any word $w$ in vocabulary $W$ of the embedding, we can define
$w_{B}$ as the part of $w$ along the gender direction.

Based on the experiments shown in Figure 10, it is justified to take the
gender direction as the (normalized) first right singular vector, $v_{B}$, or
the full data set data projected onto the subspace $F$. Then, the component of
a word vector $w$ along $v_{B}$ is simply $\langle w,v_{B}\rangle v_{B}$.

Calculating this component when the gender subspace is defined by two or more
of the top right singular vectors of $V$ can be done similarly.

We should note here that the gender subspace defined here passes through the
origin. Centering the data and using PCA to define the gender subspace lets
the gender subspace not pass through the origin. We see a comparison in the
two methods in Section 5 as $HD$ uses PCA and we use SVD to define the gender
direction.

Figure 10: Fractional singular values for (a) male-female word pairs (b) one
gendered word - one random word (c) random word pair (d) random unit vectors

##  Appendix D Word Lists

###  D.1 Word Pairs used for Flipping

actor actress
author authoress
bachelor spinster
boy girl
brave squaw
bridegroom bride
brother sister
conductor conductress
count countess
czar czarina
dad mum
daddy mummy
duke duchess
emperor empress
father mother
father-in-law mother-in-law
fiance fiancee
gentleman lady
giant giantess
god goddess
governor matron
grandfather grandmother
grandson granddaughter
he she
headmaster headmistress
heir heiress
hero heroine
him her
himself herself
host hostess
hunter huntress
husband wife
king queen
lad lass
landlord landlady
lord lady
male female
man woman
manager manageress
manservant maidservant
masseur masseuse
master mistress
mayor mayoress
milkman milkmaid
millionaire millionairess
monitor monitress
monk nun
mr mrs
murderer murderess
nephew niece
papa mama
poet poetess
policeman policewoman
postman postwoman
postmaster postmistress
priest priestess
prince princess
prophet prophetess
proprietor proprietress
shepherd shepherdess
sir madam
son daughter
son-in-law daughter-in-law
step-father step-mother
step-son step-daughter
steward stewardess
sultan sultana
tailor tailoress
uncle aunt
usher usherette
waiter waitress
washerman washerwoman
widower widow
wizard witch

###  D.2 Occupation Words

detective
ambassador
coach
officer
epidemiologist
rabbi
ballplayer
secretary
actress
manager
scientist
cardiologist
actor
industrialist
welder
biologist
undersecretary
captain
economist
politician
baron
pollster
environmentalist
photographer
mediator
character
housewife
jeweler
physicist
hitman
geologist
painter
employee
stockbroker
footballer
tycoon
dad
patrolman
chancellor
advocate
bureaucrat
strategist
pathologist
psychologist
campaigner
magistrate
judge
illustrator
surgeon
nurse
missionary
stylist
solicitor
scholar
naturalist
artist
mathematician
businesswoman
investigator
curator
soloist
servant
broadcaster
fisherman
landlord
housekeeper
crooner
archaeologist
teenager
councilman
attorney
choreographer
principal
parishioner
therapist
administrator
skipper
aide
chef
gangster
astronomer
educator
lawyer
midfielder
evangelist
novelist
senator
collector
goalkeeper
singer
acquaintance
preacher
trumpeter
colonel
trooper
understudy
paralegal
philosopher
councilor
violinist
priest
cellist
hooker
jurist
commentator
gardener
journalist
warrior
cameraman
wrestler
hairdresser
lawmaker
psychiatrist
clerk
writer
handyman
broker
boss
lieutenant
neurosurgeon
protagonist
sculptor
nanny
teacher
homemaker
cop
planner
laborer
programmer
philanthropist
waiter
barrister
trader
swimmer
adventurer
monk
bookkeeper
radiologist
columnist
banker
neurologist
barber
policeman
assassin
marshal
waitress
artiste
playwright
electrician
student
deputy
researcher
caretaker
ranger
lyricist
entrepreneur
sailor
dancer
composer
president
dean
comic
medic
legislator
salesman
observer
pundit
maid
archbishop
firefighter
vocalist
tutor
proprietor
restaurateur
editor
saint
butler
prosecutor
sergeant
realtor
commissioner
narrator
conductor
historian
citizen
worker
pastor
serviceman
filmmaker
sportswriter
poet
dentist
statesman
minister
dermatologist
technician
nun
instructor
alderman
analyst
chaplain
inventor
lifeguard
bodyguard
bartender
surveyor
consultant
athlete
cartoonist
negotiator
promoter
socialite
architect
mechanic
entertainer
counselor
janitor
firebrand
sportsman
anthropologist
performer
crusader
envoy
trucker
publicist
commander
professor
critic
comedian
receptionist
financier
valedictorian
inspector
steward
confesses
bishop
shopkeeper
ballerina
diplomat
parliamentarian
author
sociologist
photojournalist
guitarist
butcher
mobster
drummer
astronaut
protester
custodian
maestro
pianist
pharmacist
chemist
pediatrician
lecturer
foreman
cleric
musician
cabbie
fireman
farmer
headmaster
soldier
carpenter
substitute
director
cinematographer
warden
marksman
congressman
prisoner
librarian
magician
screenwriter
provost
saxophonist
plumber
correspondent
organist
baker
doctor
constable
treasurer
superintendent
boxer
physician
infielder
businessman
protege

###  D.3 Names used for Gender Bias Detection

Male : { john, william, george, liam, andrew, michael, louis, tony, scott,
jackson }
Female : { mary, victoria, carolina, maria, anne, kelly, marie, anna, sarah,
jane }

###  D.4 Names used for Racial Bias Detection and Dampening

European American : { brad, brendan, geoffrey, greg, brett, matthew, neil,
todd, nancy, amanda, emily, rachel }
African American : { darnell, hakim, jermaine, kareem, jamal, leroy, tyrone,
rasheed, yvette, malika, latonya, jasmine }
Hispanic : { alejandro, pancho, bernardo, pedro, octavio, rodrigo, ricardo,
augusto, carmen, katia, marcella , sofia }

###  D.5 Names used for Age related Bias Detection and Dampening

Aged : { ruth, william, horace, mary, susie, amy, john, henry, edward,
elizabeth }
Youth : { taylor, jamie, daniel, aubrey, alison, miranda, jacob, arthur,
aaron, ethan }


# Backward energy flow in simple 4-wave electromagnetic fields

Peeter Saari1,2 and Ioannis Besieris3 1Institute of Physics, University of
Tartu, W. Ostwaldi 1, 50411, Tartu, Estonia  2Estonian Academy of Sciences,
Kohtu 6, 10130 Tallinn, Estonia  3The Bradley Department of Electrical and
Computer Engineering, Virginia Polytechnic Institute and State University,
Blacksburg, Virginia 24060, USA  Corresponding author<EMAIL_ADDRESS>
###### Abstract

Electromagnetic energy backflow is a phenomenon occurring in regions where the
direction of the Poynting vector is opposite to that of the propagation of the
wave field. It is particularly remarkable in the nonparaxial regime and has
been exhibited in the focal region of sharply focused beams, for vector Bessel
beams, and vector-valued spatiotemporally localized waves. A detailed study is
undertaken of this phenomenon and the conditions for its appearance are
examined in detail in the case of a superposition of four plane waves in free
space, the simplest electromagnetic arrangement for the observation of
negative energy flow, as well as its comprehensive and transparent physical
interpretation. It is shown that the state of polarization of the constituent
components of the electromagnetic plane wave quartet determines whether energy
backflow takes place or not and what values the energy flow velocity assumes.
Depending on the polarization angles, the latter can assume any value from $c$
(the speed of light in vacuum) to $-c$ in certain spatiotemporal regions of
the field.

  * Keywords: Poynting vector, energy velocity, polarization, energy backflow

##  1 Introduction

Energy backflow of a free-space electromagnetic (EM) field is an unusual
effect exhibited in regions where the direction of the Poynting vector is
opposite to the propagation direction of the wave field. In this sense, the
terms ’reversed’ or ’negative’ Poynting vector are also used. The effect is
particularly remarkable in the case of non-paraxial propagation and occurs,
e.g., in the focal region of sharply focused optical beams [1]-[3], in
superpositions of four plane waves directed under equal angles $\theta$ with
respect to the propagation axis [4, 5], for vector Bessel beams [6]-[10] and
vector X-waves [11]. Note that Bessel beams are superpositions of infinitely
many monochromatic plane waves making an angle $\theta$ with respect to a
specific axis (say $z$) which corresponds to the direction of propagation of
the Bessel beam, and X-waves—representatives of the rich family of so-called
localized waves—are just the same superpositions of ultrashort pulses
[12]-[15]. Energy backflow has attracted much interest in connection with
applications in microparticle manipulation like ’tractor beams’, etc., see
reviews [16, 17].

Our purpose is to study in detail the mechanism and conditions for the
appearance of the reversed Poynting vector in a quartet of plane waves—the
simplest of EM fields for observation of the effect. In distinction from
earlier studies, we consider time-dependent EM fields instead of common cycle-
averaged ones, and make use of the energy flow velocity. Although the latter
is given by a well-known formula related to the Poynting vector and the energy
density, it exhibits some surprising properties. We hope that such a study
will result in a more comprehensible and transparent picture of the nature of
the backward energy flow effect than it can be deduced from cumbersome
formulas describing the effect in the case of the vector Bessel beams and
other more complicated fields.

In the next two sections we shall study a quartet of plane waves and show how
the polarizations of the constituents of the quartet determine whether the
energy backflow takes place or not, and what values the energy flow velocity
assumes. Section IV is devoted to studying of a general dependence of the
energy flow velocity on the polarization angles. In section V we discuss the
results, touch on the problem of physical interpretation of the Poynting
vector and some alternative definitions—noteworthy in the given context—of the
energy flow velocity. Concluding remarks are provided in section VI.

##  2 4-wave fields with regions of negative energy velocity

A quartet of plane waves symmetrically directed with respect to the
propagation axis of the whole field constitute the simplest EM field for which
energy backflow is possible111There is a claim in the literature [16] that a
superposition of two plane waves is sufficient to introduce a negative
Poynting vector. This is not wrong, but in this case the vector is not
reversed with respect to the propagation direction of the resultant field.. In
order to study the superposition of four interfering plane waves of different
polarizations, we start from a ”seed” plane wave propagating along the axis
$z$ with the wavevector $\mathbf{k}_{0}=\left(0,0,k\right)$ and polarization
vector $\mathbf{e}_{0}\left(\phi\right)=\left(\cos\phi,\sin\phi,0\right)$,
where $\phi$ is an arbitrary angle of (linear) polarization with respect to
the axis $x$. By making use of well-known rotation $\left(3\times
3\right)$-matrices $\mathbf{R}_{x}\left(\theta\right)$ and
$\mathbf{R}_{y}\left(\theta\right)$ which rotate vectors by an angle $\theta$
about the $x$\- and $y$-axis, respectively, we express the wave vectors of our
quartet of plane waves (see figure 1) as follows:

$\displaystyle\mathbf{k}_{1}$
$\displaystyle=\mathbf{R}_{y}\left(\theta\right)\mathbf{k}_{0}\;,$
$\displaystyle\mathbf{k}_{2}$
$\displaystyle=\mathbf{R}_{x}\left(-\theta\right)\mathbf{k}_{0}\;,$ (1)
$\displaystyle\mathbf{k}_{3}$
$\displaystyle=\mathbf{R}_{y}\left(-\theta\right)\mathbf{k}_{0}\;,$
$\displaystyle\mathbf{k}_{4}$
$\displaystyle=\mathbf{R}_{x}\left(\theta\right)\mathbf{k}_{0}\;.$

In this section we study the case where the ”seed” wave for all four waves is
polarized along the axis $x$, i.e., $\phi=0$. Consequently, we obtain the
following polarization vectors in the quartet:

$\displaystyle\mathbf{e}_{1}$
$\displaystyle=\mathbf{R}_{y}\left(\theta\right)\mathbf{e}_{0}\left(0\right)=\left(\cos\theta,0,-\sin\theta\right),$
$\displaystyle\mathbf{e}_{2}$
$\displaystyle=\mathbf{R}_{x}\left(-\theta\right)\mathbf{e}_{0}\left(0\right)=\left(1,0,0\right),$
(2) $\displaystyle\mathbf{e}_{3}$
$\displaystyle=\mathbf{R}_{y}\left(-\theta\right)\mathbf{e}_{0}\left(0\right)=\left(\cos\theta,0,\sin\theta\right),$
$\displaystyle\mathbf{e}_{4}$
$\displaystyle=\mathbf{R}_{x}\left(\theta\right)\mathbf{e}_{0}\left(0\right)=\left(1,0,0\right),$

which are also depicted in figure 1.

Figure 1: Wave and polarization vectors for the quartet of plane waves
according to Eqs. (1) and Eq. (2).

Our aim is to study not only spatial but also temporal dependence of energy
flow, i.e. we cannot restrict ourselves to calculation of a cycle-averaged
Poynting vector. Therefore, instead of complex field expressions we have to
work with real ones. We assume a spatio-temporal dependence of the constituent
plane waves of the form $\cos\left(\mathbf{kr}-kct\right)$, where
$k=\omega/c=2\pi/\lambda$ and $\mathbf{r}=\left(x,y,z\right)$. As a matter of
fact, if the cosine is multiplied by a pulse envelope or replaced by any other
suitable real function describing the pulse, some final results will still
hold. So we get for the electric field of the quartet

$\displaystyle\mathbf{E}\left(x,y,z,ct\right)$
$\displaystyle={\displaystyle\sum\limits_{j=1}^{4}}\cos\left(\mathbf{k}_{j}\mathbf{r}-kct\right)\mathbf{e}_{j}=$
(6) $\displaystyle=\left(\begin{array}[c]{c}2\cos(kz\cos\theta-
kct)\left[\cos\left(ky\sin\theta\right)+\cos\left(kx\sin\theta\right)\cos\theta\right]\\\
0\\\ 2\sin(kz\cos\theta-
kct)\sin\left(kx\sin\theta\right)\sin(\theta)\end{array}\right),$

where unity amplitude of each plane wave has been assumed. Eq. (6) indicates
that the field propagates along the axis $z$ due to the symmetry of the pairs
of plane waves in the quartet. Likewise for the magnetic field we get

$\displaystyle\mathbf{B}\left(x,y,z,ct\right)$
$\displaystyle={\displaystyle\sum\limits_{j=1}^{4}}\cos\left(\mathbf{k}_{j}\mathbf{r}-kct\right)\left(\mathbf{k}_{j}\times\mathbf{e}_{j}\right)=$
(10) $\displaystyle=\left(\begin{array}[c]{c}0\\\ 2\cos(kz\cos\theta-
kct)\left[\cos\left(kx\sin\theta\right)+\cos\left(ky\sin\theta\right)\cos\theta\right]\\\
2\sin(kz\cos\theta-
kct)\sin\left(ky\sin\theta\right)\sin(\theta)\end{array}\right),$

where Heaviside–Lorentz units ($\varepsilon_{0}=\mu_{0}=1$) and normalized
speed of light $c=1$ are assumed. In the given units the Poynting vector
$\mathbf{S}$ and the energy density $w$ are given by

$\displaystyle\mathbf{S}\left(x,y,z,ct\right)$ $\displaystyle=$
$\displaystyle\mathbf{E}\left(x,y,z,ct\right)\times\mathbf{B}\left(x,y,z,ct\right),$
$\displaystyle w\left(x,y,z,ct\right)$ $\displaystyle=$
$\displaystyle\frac{1}{2}\left[\mathbf{E}^{2}\left(x,y,z,ct\right)+\mathbf{B}^{2}\left(x,y,z,ct\right)\right].$
(11)

For the quartet the expression of $\mathbf{S}$ turns out to be not too
cumbersome if one introduces the following rescaled and dimensionless
coordinates: $\widetilde{x}\equiv k~{}x\sin\theta$, $\widetilde{y}\equiv
k~{}y\sin\theta$, $\widetilde{z}\equiv k~{}z\cos\theta$, $\widetilde{t}\equiv
kct$. With these definitions we obtain

$\mathbf{S}\left(\widetilde{x},\widetilde{y},\widetilde{z},\widetilde{t}\right)=\left(\begin{array}[c]{c}-4\cos\left(\widetilde{z}-\widetilde{t}\right)\sin\left(\widetilde{z}-\widetilde{t}\right)\sin\theta\sin\widetilde{x}~{}\left(\cos\widetilde{x}+\cos\widetilde{y}\cos\theta\right)\\\
-4\cos\left(\widetilde{z}-\widetilde{t}\right)\sin\left(\widetilde{z}-\widetilde{t}\right)\sin\theta\sin\widetilde{y}~{}\left(\cos\widetilde{y}+\cos\widetilde{x}\cos\theta\right)\\\
4\cos^{2}\left(\widetilde{z}-\widetilde{t}\right)\left(\cos\widetilde{x}+\cos\widetilde{y}\cos\theta\right)\left(\cos\widetilde{y}+\cos\widetilde{x}\cos\theta\right)\end{array}\right).$
(12)

From Eq. (12) the following conclusions can be drawn:

  1. 1.

the field of the Poynting vector propagates as a whole with velocity
$v=c/\cos\theta$ along the axis $z$ and it is invariant because it depends on
time only through the difference $z\cos\theta-ct$ (we name it as the
propagation variable);

  2. 2.

if $z=t=0$ and generally on planes $k(z\cos\theta-ct)=n\pi,~{}n\in\mathbb{Z}$
the transverse components of $\mathbf{S}$ vanish;

  3. 3.

on planes $k(z\cos\theta-ct)=(2n+1)\pi/2,~{}n\in\mathbb{Z}$ the field
$\mathbf{S}$ vanishes;

  4. 4.

the longitudinal component $S_{z}$ is invariant with respect to interchange of
coordinates $x$ and $y$;

  5. 5.

the transverse components transform to each other if the coordinates $x$ and
$y$ are interchanged;

  6. 6.

there are regions in any transverse plane where the longitudinal component
$S_{z}$ is negative, i.e., the energy flows backward and this phenomenon takes
place irrespectively of time instant or cycle averaging, although the latter
affects the numerical value of negative $S_{z}$.

The energy density is given as follows in terms of the dimensionless
coordinates:

$\displaystyle
w\left(\widetilde{x},\widetilde{y},\widetilde{z},\widetilde{t}\right)$
$\displaystyle=$ $\displaystyle
2\cos^{2}\left(\widetilde{z}-\widetilde{t}\right)[\left(1+\cos^{2}\theta\right)\cos^{2}\widetilde{x}+4\cos\theta\cos\widetilde{x}\cos\widetilde{y}+$
$\displaystyle\left(1+\cos^{2}\theta\right)\cos^{2}\widetilde{y}]+2\sin^{2}\left(\widetilde{z}-\widetilde{t}\right)\left(\sin^{2}\widetilde{x}+\sin^{2}\widetilde{y}\right)\sin^{2}\theta.$

Similarly to $S_{z}$, the energy density obeys the symmetry
$w\left(x,y,z,ct\right)=w\left(y,x,z,ct\right).$ This conclusion is confirmed
by numerically calculated plots, an example of which is presented in figure 2.

Figure 2: Energy density and energy flow for the quartet of plane EM waves
with polarization vectors given by Eq. (2) and $\theta=2\pi/5=72^{\circ}$.
(a)—greyscale plot of $w$ as a function of the transverse coordinate
$\widetilde{x}$ at $\widetilde{y}=0$ (or $\widetilde{y}$ at $\widetilde{x}=0$)
and the propagation variable $\widetilde{z}-\widetilde{t}$; (b)—the same with
decreased contrast and superimposed by a vector field plot of the energy flow
velocity. Shown are projections of the velocity vectors onto the plane
$\left(z,x\right)$; (c)—dependence on $\widetilde{x}$ (along vertical axis) of
normalized longitudinal ($z$-) component of the Poynting vector
$S_{z}(\widetilde{x},0,0,0)/S_{z}(0,0,0,0)$ showing the minimal values
$-0.28$. The range of the vertical axes is $[-2.5\pi,2.5\pi]$ and of the
propagation axes is $[-1.25\pi,1.25\pi]$. In plot (c) the scale of the
transversal coordinate $\widetilde{x}=kx\sin\theta$ along the vertical is
shown in units of $\pi$.

The interference pattern in the energy density (figure 2a) remarkably differs
from that of two plane waves: if the pair (waves # 2,4) were removed from the
quartet, all maxima would become of equal intensity (e.g., five equally bright
maxima along the central vertical $\widetilde{z}-\widetilde{t}=0$ instead of
three in figure 2a). Since the magnitude of the Poynting vector $\mathbf{S}$
is strong enough only in regions of local maxima and therefore a vector field
plot of $\mathbf{S}$ would be of low distinctness, we present the vector field
plot of the energy flow velocity $\mathbf{V}=\mathbf{S}/w$ instead. Figure 2b
clearly reveals that at locations of weak maxima of $\ w$ the energy flows
backward. While the strength of the component $S_{z}$ at these locations
comprises only 28 % (see figure 2c) of that in locations of strong maxima of
$w$, the velocities $\mathbf{V}$ have exactly the same absolute value 1 (or
$c$ with dimensions) but the opposite signs. Numerical evaluation shows—and it
is evident from figure 2b also—that wherever the velocity is purely
longitudinal it equals to the physically maximal possible value $\pm c$ (small
irregularities in the vector field at points of vanishing $w$ are
computational artifacts due to ”zero divide by zero”). We point out the very
noteworthy result that energy backflow occurs also with velocity $c$. While
according to Eqs. (6) and (10) neither of the fields $\mathbf{E}$ or
$\mathbf{B}$ is transverse, at locations where $V_{z}\mathbf{=}\pm c$ the
fields are locally TEM, mutually perpendicular and of equal magnitude (in the
chosen units)—as they have to be for a null EM field.

Note that the velocity components $V_{z}$ and $V_{x}$ on the plane
$\left(y=0\right)$ constitute the vector field plot in figure 2b. If instead
we plotted figure 2b with $V_{z}$ and $V_{y}$ in the same plane, we would
obtain another projection of the energy flow velocity, where all arrows would
be directed with respect to the axis $z$ either parallel (in the regions of
stronger maxima of $w$) or antiparallel (in the regions of weaker maxima of
$w$). This is understandable because according to Eq. (12) $S_{y}$ vanishes on
the plane $\left(y=0\right)$.

Naturally, the smaller the angle $\theta$ the less remarkable will be the
backflow effect: peak negative values of the normalized $S_{z}$ decrease to
$11\%$ if $\theta=\pi/3$ and to $3\%$ if $\theta=\pi/4$. Of course, in the
case of paraxial geometry the energy backflow vanishes. However, the
abovementioned characteristic features of the velocity vector field remain
unchanged with decreasing $\theta$, except transverse narrowing of the
backflow channels.

Figure 3: Dependence of the longitudinal component of the energy velocity
$V_{z}$ (solid curve) and the normalized energy density (dashed curve)
$\widetilde{w}\equiv
w(\widetilde{x},\widetilde{y},\widetilde{z},\widetilde{t})/w(0,0,0,0)$ on the
dimensionless transverse coordinate $\widetilde{x}\equiv kx\sin\theta$;
(a)—along the line $y=0$ on the plane $(\widetilde{z}-\widetilde{t}=0)$ and
(b)—along the line $y=0$ on the plane $(\widetilde{z}-\widetilde{t}=0.4\pi)$.
Note that in the plot (b) the scale for $\widetilde{w}$ is on the right-hand
side and differs from the joint scale for $V_{z}$ and $\widetilde{w}$ in (a).

Similar narrowing of the backflow channels is observable also if one moves
away from the planes $\widetilde{z}-\widetilde{t}=n\pi,~{}n\in\mathbb{Z}$.
Curves in figure 3 show the dependence of the longitudinal component of the
energy flow velocity in comparison with that of the energy density along the
central vertical in figure 2b, i.e., where $\widetilde{z}-\widetilde{t}=0$. We
see in figure 3a that transverse widths of the regions where $V_{z}=\pm c$ are
approximately equal to FWHM of the maxima of $w$. But on the plane
$\widetilde{z}-\widetilde{t}=0.4\pi$, which is close to the plane
$\widetilde{z}-\widetilde{t}=\pi/2$ where $V=0$, the backflow regions narrow
substantially and their FWHM become about two times less than those of the
minima or maxima of the energy density. By introducing the separate secondary
scale for
$\widetilde{w}\left(x,0,0.4\pi,0\right)=w\left(x,0,0.4\pi,0\right)/w\left(0,0,0,0\right)$
on the right border of the plot, it is clearly seen in figure 3b that the
negative minima of $V_{z}$ are encompassed by the (positive) minima of
$\widetilde{w}$.

In the discussion so far it has been assumed that the seed polarization angle
$\phi=0$. However, it turns out in general that neither $\mathbf{S}$ nor
$\mathbf{V}$ depend on the angle of polarization $\phi$ (if its value is the
same for all four waves). This is understandable because despite the fact that
the polarization vectors of the quartet do depend on $\phi$ according to Eq.
(2), the combination of their mutual relative directions, which determines
$\mathbf{S}$ and $w$, is invariant.

##  3 4-wave fields without regions of negative energy velocity

A closer study of the effect of the angle $\phi$ on the energy flow reveals
that what is essential is the difference between the angle $\phi_{13}$ for the
pair $(1,3)$ and the angle $\phi_{24}$ for the pair $(2,4)$ in the quartet of
waves. Thus, if $\phi_{13}-\phi_{24}=0$, we return to the results of the
previous section. Of special interest are the cases
$\phi_{13}-\phi_{24}=\pm\pi/2$. Since the polarization vectors together with
fields $\mathbf{E}$ and $\mathbf{B}$ depend on both $\phi_{13}$ and
$\phi_{24}$, we prescribe values to them separately and for simplicity make
one of them equal to zero. So, first we consider the case $\phi_{13}=\pi/2$
and $\phi_{24}=0$. Then $\mathbf{e}_{1}=\mathbf{e}_{3}=\left(0,1,0\right)$ and
$\mathbf{e}_{2}=\mathbf{e}_{4}=\left(1,0,0\right),$ i.e., the pair $(1,3)$ is
nothing but a replica of the pair $(2,4)$ rotated about the axis $z$ clockwise
by $\pi/2$. Analogously to Eq. (6) and introducing the dimensionless
arguments, we obtain

$\mathbf{E}\left(\widetilde{x},\widetilde{y},\widetilde{z},\widetilde{t}\right)=\left(\begin{array}[c]{c}2\cos\left(\widetilde{z}-\widetilde{t}\right)\cos\widetilde{y}\\\
2\cos\left(\widetilde{z}-\widetilde{t}\right)\cos\widetilde{x}\\\
0\end{array}\right)$ (13)

and likewise for the magnetic field

$\mathbf{B}\left(\widetilde{x},\widetilde{y},\widetilde{z},\widetilde{t}\right)=\left(\begin{array}[c]{c}-2\cos\left(\widetilde{z}-\widetilde{t}\right)\cos\widetilde{x}\cos\theta\\\
2\cos\left(\widetilde{z}-\widetilde{t}\right)\cos\widetilde{y}\cos\theta\\\
-2\sin\left(\widetilde{z}-\widetilde{t}\right)\left(\sin\widetilde{x}-\sin\widetilde{y}\right)\sin\theta\end{array}\right).$
(14)

We see that in the given case the EM field as a whole is transverse electric
(TE). Eq. (11) gives for the Poynting vector and the energy density the
expressions

$\mathbf{S}\left(\widetilde{x},\widetilde{y},\widetilde{z},\widetilde{t}\right)=\left(\begin{array}[c]{c}-4\cos\left(\widetilde{z}-\widetilde{t}\right)\sin\left(\widetilde{z}-\widetilde{t}\right)\sin\theta\cos\widetilde{x}~{}\left(\sin\widetilde{x}-\sin\widetilde{y}\right)\\\
4\cos\left(\widetilde{z}-\widetilde{t}\right)\sin\left(\widetilde{z}-\widetilde{t}\right)\sin\theta\cos\widetilde{y}~{}\left(\sin\widetilde{x}-\sin\widetilde{y}\right)\\\
4\cos^{2}\left(\widetilde{z}-\widetilde{t}\right)\cos\theta\left(\cos^{2}\widetilde{x}+\cos^{2}\widetilde{y}\right)\end{array}\right)$
(15)

and

$\displaystyle
w\left(\widetilde{x},\widetilde{y},\widetilde{z},\widetilde{t}\right)$
$\displaystyle=$ $\displaystyle
2\cos^{2}\left(\widetilde{z}-\widetilde{t}\right)\left(1+\cos^{2}\theta\right)\left(\cos^{2}\widetilde{x}+\cos^{2}\widetilde{y}\right)+$
$\displaystyle
2\sin^{2}\left(\widetilde{z}-\widetilde{t}\right)\sin^{2}\theta\left(\sin^{2}\widetilde{x}-\sin^{2}\widetilde{y}\right).$

Conclusions from Eq. (15) are the same as those derived from Eq. (12) in the
previous section, but with one essential exception: the 6th conclusion now is
that the longitudinal component $S_{z}$ cannot take negative values, i.e., the
energy flows only forward (of course, it is assumed that $|\theta|\leq\pi/2$
since otherwise the constituent plane waves would propagate backward). These
conclusions, including the absence of energy backward flow are confirmed by
numerically calculated plots, see figure 4.

Figure 4: Energy density (a) and energy flow (b,c) of the quartet of EM plane
waves with polarization vectors
$\mathbf{e}_{1}=\mathbf{e}_{3}=\left(0,1,0\right)$ and
$\mathbf{e}_{2}=\mathbf{e}_{4}=\left(1,0,0\right)$ in which case energy flows
only in the positive direction of the axis $z$. For other parameters and
scales see the caption of figure 2.

As seen in figure 4, all maxima of the energy density are of equal intensity,
the energy velocity vector field is the same at all maxima and contains no
backward flow. But what is very noteworthy is that everywhere the energy
velocity is essentially subluminal (less than $c$) and wherever the energy
flows along the propagation axis $z,$ its velocity is given by the expression

$\frac{\left|\mathbf{V}\right|}{c}=\frac{V_{z}}{c}=\frac{2\beta}{\beta^{2}+1}=\frac{2\beta^{-1}}{\beta^{-2}+1}~{},$
(16)

where $\beta=\cos\theta$. In figure 4, we obtain for $\theta=2\pi/5$ the value
$V_{z}=0.564c$. Precisely the same value comes out from Eq. (16). For
comparison, as mentioned in the previous section, all fields $\mathbf{E,B,S}$,
and $w$ propagate with velocity $v=c/\cos\theta=\beta^{-1}c=3.24c,$ i.e.,
superluminally. Eq. (16) was first found in [18] where it was shown that
various fields—crossing plane waves, Bessel beams, Bessel-X pulses and other
propagation-invariant wave fields obey this relation between the axial phase-
or group velocity and energy flow velocity.

Next, we consider the case $\phi_{13}=0$ and $\phi_{24}=\pi/2$. Then the
polarization vectors in the quartet turn out to be

$\displaystyle\mathbf{e}_{1}$
$\displaystyle=\left(\cos\theta,0,-\sin\theta\right),$
$\displaystyle\mathbf{e}_{2}$
$\displaystyle=\left(0,\cos\theta,-\sin\theta\right),$
$\displaystyle\mathbf{e}_{3}$
$\displaystyle=\left(\cos\theta,0,\sin\theta\right),$ (17)
$\displaystyle\mathbf{e}_{4}$
$\displaystyle=\left(0,\cos\theta,\sin\theta\right),$

i.e., the pair $(2,4)$ is nothing but a replica of the pair $(1,3)$ rotated
about the axis $z$ clockwise by $\pi/2$. In contradistinction to the previous
case, the EM field now is transverse magnetic (TM) and the magnetic field
turns out to be equal to Eq. (13) but with minus sign in front of the first
component. The expression for the Poynting vector turns out to be the same as
in Eq. (15) if in the latter the difference of sines is replaced by their sum
and the factor $4$ in the second component by $-4$. As a result, the quartet
with polarization vectors of Eq. (17) exhibits the same subluminal velocity of
the energy flow as in the previous case and figure 4 illustrates equally well
the present case. In [18] we have shown that (i) a TM field of a symmetrical
pair of plane waves—the propagating direction of the first one lies on the
$(x,z)$ plane and is inclined by angle $+\theta$ with respect to the $z$-axis,
and the second one by angle $-\theta$ on the same plane—exhibit subluminal
energy velocity along the axis $z$ according to Eq. (16); (ii) any
superpositions of such pairs rotated about the axis $z$ have the same
property.. The quartet of the present case is just a particular example of
such superposition where the two pairs are rotated clockwise by $\pi/2$.

To conclude, from the results of sections 2 and 3 a general rule follows: if
the EM field of the wave quartet is purely TE or TM, there is no backward
energy flow and its (forward) velocity is subluminal obeying Eq. (12);
otherwise not only forward but also backward energy flow takes place, both
with luminal velocity $c$.

##  4 Dependence of energy velocity on polarization angles: general case

Since the energy flow velocity field so drastically depends on whether
$\phi_{13}=\pi/2$ or $\phi_{13}=0$, it is of interest to study the dependence
over the full range of angle $\phi_{13}=\phi\in[0,\pi]$ (keeping $\phi_{24}\
=0)$. Let us focus, in particular, on regions where $\mathbf{V}$ is parallel
or antiparallel to the axis $z$ and has magnitude
$\left|\mathbf{V}\right|/c=1$ or $V_{z}/c=\pm 1$, see figure 2. For the sake
of definiteness we consider two points on the axis $x$: first, the point where
$\widetilde{x}=\widetilde{y}=\widetilde{z}-\widetilde{t}=0$ and, second, where
$\widetilde{x}=\pi$, $\widetilde{y}=\widetilde{z}-\widetilde{t}=0$. For the
first point we obtain

$\frac{V_{z}}{c}=\frac{2\cos\theta+\cos\phi+\cos\phi\cos^{2}\theta}{\cos^{2}\theta+1+2\cos\phi\cos\theta},$
(18)

which for $\phi=0$ reduces to $V_{z}/c=1$ as it has to, see figure 2 and for
$\phi=\pi/2$ reduces to Eq. (16) as it also has to, see figure 4. For the
second point a slightly modified version of Eq. (18) comes out: only a
replacement $\cos\phi\rightarrow-\cos\phi$ must be done. Again, the modified
expression for $\phi=0$ reduces to $V_{z}/c=-1$ and for $\phi=\pi/2$ to Eq.
(16), as expected. The behaviour of $V_{z}/c$ over the full range of angle is
depicted in figure 5.

Figure 5: Dependence of the energy flow velocity on the angle $\phi$ for
points $\widetilde{x}=0$ (solid curve) and $\widetilde{x}=\pi$ (dashed curve)
according to Eq. (18).

We see that if $\phi=\pi/2$ the velocity is given by Eq. (16) in accord with
figure 4. At values $\phi$ close to $\pi$ the regions of forward and backward
flow of energy interchange. Equating the nominator of Eq. (18) to zero allows
one to answer the question: what are the values of $\phi$ at which the energy
backflow disappears, i.e., both curves in figure 5 become positive? The answer
is: $\phi$ must obey the condition $\cos\phi=\pm
2\cos\theta/\left(\cos^{2}\theta+1\right)$, i.e., the right-hand-side of Eq.
(16). Whether such a surprising coincidence has a deeper physical meaning,
needs further study.

Similarly to the preceding section we also studied the behaviour of the
velocity if the values of angles $\phi_{13}$ and $\phi_{24}$ are interchanged,
i.e. the velocity dependence on $\phi_{24}$ over the full range of angle
$\phi_{24}=\phi\in[0,\pi]$ (keeping $\phi_{13}\ =0)$. Although the interchange
of the angles alters the fields $\mathbf{E}$ and $\mathbf{B}$ (interchange
their components and signs), the Poynting vector $\mathbf{S}$, as well as the
energy density $w$ and therefore also $\mathbf{V}$ remain unchanged. Hence,
the results given above hold in this case as well.

##  5 Discussion

Since the variable along the horizontal axis in figures 2a, 2b, 4a, and 4b is
$z\cos\theta-ct$, these plots can be viewed either as $z$-dependence or time
dependence of the fields. The latter viewpoint allows one to read out from
these plots what would be the temporally cycle-averaged quantities commonly
evaluated from complex-valued EM fields. We see that the direction of energy
flow does not flip in time (and/or along the axis $z$) and, in particular,
where it flows backward in figure 2b, it does so all time. Also, the
transverse components of the velocity vanish as a result of such averaging.
This follows already from symmetry considerations and from the observation
that the product
$\cos\left(\widetilde{z}-\widetilde{t}\right)\sin\left(\widetilde{z}-\widetilde{t}\right)$
entering the transverse components in Eqs. (12) and (15) vanishes after cycle
averaging. In Ref. [5] the cycle-averaged $z$-component $P_{z}$ of the
Poynting vector has been evaluated for the same polarization vectors as in our
Eq. (2) and plotted on the plane $(z=0)$. Our more detailed results are in
accord with corresponding results in Ref. [5]. An example of a wave quartet
which does not reveal negative $P_{z}$ has been also studied by the authors of
the cited article, but they did it with another geometry of the polarization
vectors than in our two cases of section 3. In their quartet one pair of the
polarization vectors is also a replica of the other rotated about the axis $z$
clockwise by $\pi/2$; the field is also TM and the cycle-averaged
$z$-component of the Poynting vector is positive everywhere. We studied the
cycle non-averaged vector field of energy velocity in such quartet and found
that in contradistinction to figure 4b, the velocity vanishes at the maxima of
the energy density $w$, whereas its longitudinal component at the minima of
$w$ is given by Eq. (16).

The general conclusion from the preceding section is in accord with what has
been found for the Bessel beams and X-waves [7, 11]: without superposition of
TE and TM fields, there is no reversed Poynting vector. This is
understandable, since the quartet is an elementary example of superpositions
that construct these much more complicated waves.

Although we considered here only monochromatic waves, the results obtained
hold also for ultrashort propagation-invariant pulses, more exactly—for apex
regions of their X-like or double conical field profile. This can be readily
proved by applying the approach used in [18, 19]] for the evaluation of the
Poynting vector of few-cycle and single-cycle light sheets. Keeping in mind
that the velocity $v=c/\cos\theta$ for such pulses is simultaneously phase and
group velocity, we encounter a striking conclusion: according to Eq. (16) the
energy flows subluminally while the pulse itself propagates superluminally!
Definitely the statement “if an energy density is associated with the
magnitude of the wave… the transport of energy occurs with the group velocity,
since that is the rate of which the pulse travels along” (Ref. [20], sec. 7.8)
cannot hold if the group velocity exceeds $c$. A possibility to resolve this
contradiction is offered by the notion of reactive energy [19].

It is well known that the Poynting vector is not defined uniquely by the
Poynting theorem. Could it be that, consequently, the energy-flow velocity we
have dealt with throughout this paper is also not defined uniquely and this is
the reason of the energy backflow effect and non-equality of its velocity to
the group velocity? The discussion of whether the Poynting vector gives the
local energy flow or not has a century-long history, see, e.g., a review [21].
In particular, an alternative definition of the Poynting vector has been
proposed which is more in harmony with the notion of group velocity [22].
However, we are not going to dig into the problem here and simply refer to
[20] (sec. 6.7), where it is stated that if one takes into account that the
Poynting vector defines also the momentum density of the field and the
definition must not violate the theory of relativity, the common expression
for the Poynting vector is unique.

An intriguing alternative definition of the energy velocity is proposed in
[23]. It springs forth from an idea that a correct velocity in a given
location must be equal to the velocity of an inertial frame in which energy
flux turns to be zero at the given location. Such approach leads to the
following relation between the magnitude
$V\equiv\left|\mathbf{V}\right|=\left|\mathbf{S}\right|/w$ of the commonly
defined energy velocity and its alternative version $V_{a}$

$\frac{V}{c}=\frac{2\left(V_{a}/c\right)}{1+\left(V_{a}/c\right)^{2}}~{}.$
(19)

This relation, which is Eq. (11) from Ref. [23] rewritten in our designations,
surprisingly coincides with Eq. (16) if we equate $V_{a}/c$ to the normalized
group velocity $\beta^{-1}$ or to its reciprocal. Indeed, in the reference
frame moving with $V_{a}=\beta$ along the axis $z$ the angle $\theta$
transforms to $\pi/2$ not only for the wave quartet but also for Bessel beams
and X-type pulsed localized waves [24] resulting according to Eq. (15) in
vanishing of $S_{z}$. However, further analysis of the alternative definition
of the energy flow velocity is outside the scope of this study.

##  6 Conclusion

The results of this study of the electromagnetic field of four interfering
plane waves demonstrate in detail how crucial is the polarization of the waves
for the existence of regions where the energy flows backward with respect to
the propagation direction of the field itself. Moreover, by adjusting the
relative angle between linear polarizations of two pairs in the quartet of
waves, one can obtain any value between $-c$ and $c$ of the energy flow
velocity at maxima or minima of energy density.

Generalizations of the results for interfering few-cycle ultrashort pulses is
straightforward. They may also be useful in studies of quantum backflow and
optical analogues of this intriguing effect in particle physics [25]-[27].

## References

## References

  * [1] Richards B and Wolf E 1959 Electromagnetic diffraction in optical systems, II. Structure of the image field in an aplanatic system Proc. R. Soc. A 253 358-79
  * [2] Kotlyar V V, Stafeev S S, Nalimov A G, Kovalev A A, Porfirev A P 2020 Mechanism of formation of an inverse energy flow in a sharp focus Phys. Rev. A 101 033811
  * [3] Li H, Wang C, Tang, M and Xinzhong A 2020 Controlled negative energy flow in the focus of a radial polarized optical beam Opt. Express 28 18607
  * [4] Katsenelenbaum B Z 1997 What is the direction of the Poynting vector? J. Commun. Technol. Electron. 42 119
  * [5] You X-L, Li C-F 2020 Dependence of Poynting vector on state of polarization arXiv:2009.04119
  * [6] Turunen J, Friberg A T 1993 Self-imaging and propagation-invariance in electromagnetic fields Pure Appl. Opt. 2 51
  * [7] Novitsky A V, Novitsky D V 2007 Negative propagation of vector Bessel beams J. Opt. Soc. Am. A 24 2844
  * [8] Sukhov S, Dogariu A 2010 On the concept of tractor beams Opt. Lett. 35 3847
  * [9] Novitsky A V, Qiu C W, Wang H 2011 Single Gradientless Light Beam Drags Particles as Tractor Beams Phys. Rev. Lett. 107 203601
  * [10] Chen J, Jack Ng J, Lin Z and Chan C T 2011 Optical pulling force Nat. Photonics 5 531
  * [11] Salem M A and Bağcı H 2011 Energy flow characteristics of vector X-Waves Opt. Express 19, 8526
  * [12] Non-Diffracting Waves 2013 edited by Hernandez-Figueroa H E, Recami E and Zamboni-Rached M (New York: J. Wiley)
  * [13] Besieris I, Abdel-Rahman M, Shaarawi A, and Chatzipetros A 1998 Two fundamental representa-tions of localized pulse solutions to the scalar wave equation Prog. Electromagn. Res.19 1
  * [14] Saari P and Reivelt 2004 Generation and classification of localized waves by Lorentz transfor-mations in Fourier space Phys. Rev. E 69 036612
  * [15] Kondakci H E and Abouraddy A F 2017 Diffraction-free space-time light sheets Nat. Phot. 11 733
  * [16] Qiu C W, Palima D, Novitsky A, Gao D, Ding W, Zhukovsky S V, Gluckstad J 2014 Engineering light-matter interaction for emerging optical manipulation applications Nanophotonics 3 181
  * [17] Li H, Cao Y, Zhou L-M, Xu X, Zhu T, Shi Y, Qiu C-W, Ding W 2020 Optical pulling forces and their applications Adv. Opt. Photon. 2 288
  * [18] Saari P, Rebane O and Besieris I 2019 Energy-flow velocities of nondiffracting localized waves Phys. Rev. A 100 013849
  * [19] Saari P and Besieris I M 2020 Reactive energy in nondiffracting localized waves Phys. Rev. A 101 023812
  * [20] Jackson J D 1998 Classical Electrodynamics, 3rd ed. (Wiley, New York)
  * [21] Gough W, 1982 Poynting in the wrong direction? Eur. J. Phys. 3 83
  * [22] Hines C O 1952 Electromagnetic energy density and flux Canadian J. of Phys. 30 123
  * [23] Johns O D 2020 Relativistically correct electromagnetic energy flow arXiv:2010.10263
  * [24] Saari P and Besieris I M 2020 Relativistic aberration and null Doppler shift within the framework of superluminal and subluminal nondiffracting waves J. Phys. Commun. 4 105011
  * [25] Berry M V 2010 Quantum backflow, negative kinetic energy, and optical retro-propagation J. Phys. A 43 415302
  * [26] Villanueva A A D 2020 The negative flow of probability Am. J. Phys 88 325
  * [27] Eliezer Y, Zacharias T and Bahabad A 2020 Observation of optical backflow Optica 7 72