# Low Row Rank Parity Check Codes

Ermes Franch  Department of Informatics, University of Bergen, Norway  Chunlei
Li  Department of Informatics, University of Bergen, Norway

(January 2021)

###### Abstract

In this paper we present an extended variant of low rank parity check matrix
(LRPC) codes that have received significant interests in recent years. It is
shown that the extension indeed yields a super-family of LRPC codes, which are
termed low row rank parity check codes. The decoding method of the proposed
codes is also investigated.

##  1 Introduction

Rank metric codes were independently introduced by Delsarte [1], Gabidulin [2]
and Roth [3]. The Gabidulin codes, which are an equivalent of the Reed-Solomon
codes for the rank metric, have found applications in network coding [4],
space-time coding [5] and cryptography [6]. In 2013 Gaborit, Murat and Zémor
introduced a new family of rank metric codes, termed low rank parity check
(LRPC) codes, which can be seen as the equivalent of LDPC codes [7]. With an
efficient probabilistic decoding algorithm, LRPC codes have been applied in
several cryptosystems [8, 9, 10] as well as a signature scheme [11], which are
designed for better resistance against quantum computers and smaller key size
compared to the classic McEliece cryptosystem in Hamming metric [12]. Among
these schemes, the cryptosystem ROLLO [10], based on a subfamily of LRPC
codes, has reached the second round of the National Institute of Standards and
Technology (NIST) post-quantum competition.

In this paper we will present an extended variant of the LRPC codes, which are
called low row rank parity check codes, by applying the constraint of low rank
weight only on each row of the parity check matrix instead of on the whole
parity check matrix as for LRPC codes. The proposed low rank parity check
codes are shown to contain LRPC codes as a strictly proper subclass. Due to
the loosened constraint, the decoding technique of LRPC codes no longer works
for the low row rank parity check codes. We investigate the decoding process
of these new codes and discuss potential reasons for decoding failures.

##  2 Background on rank metric codes

###  2.1 Notations

Let $\mathbb{F}_{q}$ the finite field of $q$ elements, where $q$ is a power of
a prime $p$. The extension field $\mathbb{F}_{q^{m}}$ can be deemed as an
$\mathbb{F}_{q}$-vector subspace of dimension $m$. Given elements
$v_{1},\ldots,v_{r}$ in $\mathbb{F}_{q^{m}}$, the $\mathbb{F}_{q}$-vector
subspace generated by $v_{i}$’s, denoted in the calligraphic upper case, is
given by

$\mathcal{V}=\langle
v_{1},\ldots,v_{r}\rangle=\left\\{\sum_{i=1}^{r}\lambda_{i}v_{i}\,|\,\lambda_{i}\in\mathbb{F}_{q}\right\\}.$

We will denote vectors over $\mathbb{F}_{q^{m}}$ in bold lower case as
$\mathbf{v}=(v_{1},\ldots,v_{n})$ and denote matrices in upper case as
$M=(M_{i,j})$. Given a $k\times n$ matrix $M$ and two non-empty ordered
subsets $I\subseteq\\{1,\ldots,k\\},J\subseteq\\{1,\ldots,n\\}$, the sub-
matrix of $M$ whose rows are indexed by $I$ and columns are indexed by $J$ is
denoted as $M_{I,J}$. The transpose of the matrix $M$ will be denoted as
$M^{\intercal}$.

Let $\mathcal{V}=\langle v_{1},\ldots,v_{d}\rangle$ be an
$\mathbb{F}_{q}$-vector space of dimension $d$. We denote by $\mathcal{V}^{*}$
the set $\mathcal{V}\setminus\\{0\\}$ and by $\mathcal{V}^{-1}$ the set
$\\{v^{-1}|v\in\mathcal{V}^{*}\\}$. Given another $\mathbb{F}_{q}$-vector
space $\mathcal{U}=\langle u_{1},\ldots,u_{d}\rangle$ of dimension $r$, the
product vector space of $\mathcal{V}$ and $\mathcal{U}$ is defined as
$\mathcal{V}.\mathcal{U}=\langle
v_{1}u_{1},\ldots,v_{1}u_{d},\ldots,v_{r}u_{1},\ldots,v_{r}u_{d}\rangle,$
which has dimension at most $rd$.

###  2.2 Basic definitions

We start with the definitions of rank weight and rank distance of elements in
$\mathbb{F}_{q^{m}}$.

######  Definition 1.

Let $S\subseteq\mathbb{F}_{q^{m}}$ and let $\langle S\rangle$ be the
$\mathbb{F}_{q}$-vector subspace generated by all the elements of $S$. The
rank weight of $S$ is given by $\mathrm{rw}(S)=\dim(\langle S\rangle)$.

We can extend this definition in a natural way to vectors and matrices.

######  Definition 2.

Let $\mathbf{v}=(v_{1},\ldots,v_{n})$ be a vector over $\mathbb{F}_{q^{m}}$.
The $\mathbb{F}_{q}$-vector subspace $\mathcal{V}=\langle
v_{1},\ldots,v_{n}\rangle$ is called the support of $\mathbf{v}$. The rank
weight of $\mathbf{v}$ is given by

$\mathrm{rw}(\mathbf{v}):=\dim\left(\mathcal{V}\right)=\dim(\langle
v_{1},\ldots,v_{n}\rangle).$

The rank distance between two vectors $\mathbf{u},\mathbf{v}$ is given by

$\mathrm{d}(\mathbf{u},\mathbf{v}):=\mathrm{rw}(\mathbf{u}-\mathbf{v}).$

######  Definition 3.

The weight $\delta$ of a matrix $H$ over $\mathbb{F}_{q^{m}}$ is defined as
the rank weight of the set of all its elements. The matrix $H$ is said to have
low weight if $\delta<m$, particularly $\delta\ll m$. A rank metric code
$\mathcal{C}$ over $\mathbb{F}_{q^{m}}$ is said to be an LRPC code of weight
$\delta$ if $\mathcal{C}$ admits a parity check matrix $H$ of weight less than
or equal to $\delta$.

Given an LRPC code $\mathcal{C}$ with a low-weight parity check matrix $H$, it
is easy to build another parity check matrix $H^{\prime}$ for $\mathcal{C}$ of
full weight $m$. On the contrary, starting from a matrix $H^{\prime}$ of
weight $m$, it is hard to find an equivalent matrix $H$ of low weight. It
would be at least as hard as finding $n-k$ linearly independent vectors of low
weight in the dual code of $\mathcal{C}$ (which are not guaranteed to exist if
$\mathcal{C}$ is not an LRPC code).

##  3 The row-LRPC codes

In this section we will introduce a super-family of LRPC codes by applying the
low-weight constraint only to rows of a parity check matrix instead of the
whole matrix. The decoding process of those new codes will be examined in the
next section.

######  Definition 4.

def:row-LRPC The row weight $\rho$ of an $(n-k)\times n$ matrix $H$ over
$\mathbb{F}_{q^{m}}$ is defined as the maximum rank weight of all its rows:

$\rho=\max\limits_{1\leq i\leq n-k}(\mathrm{rw}(\mathbf{h}_{i})),$

where $\mathbf{h}_{i}$ for $1\leq i\leq n-k$ is the $i$-th row of $H$. The
matrix $H$ is said to have low row weight if $\rho<m$, particularly $\rho\ll
m$. We call a rank metric code $\mathcal{C}$ a low row rank parity check
(shortly row-LRPC) code of row weight $\rho$ if $\mathcal{C}$ admits a parity
check matrix $H$ of row weight less than or equal to $\rho$.

In the above definition we use the abbreviation row-LRPC codes, instead of
LRRPC codes, for better distinguishing them from LRPC codes. For simplicity of
presentation, we denote by $\mathrm{LRPC}_{q}(n,m,k,\delta)$ the set of all
the $[n,k]$ linear LRPC codes over $\mathbb{F}_{q^{m}}$ of weight $\delta$,
and by row-$\mathrm{LRPC}_{q}(n,m,k,\rho)$ the set of all the $[n,k]$ row-LRPC
codes over $\mathbb{F}_{q^{m}}$ of row weight $\rho$.

######  Remark 1.

From Definition 3 and Definition LABEL:def:row-LRPC, we have
$\mathrm{LRPC}_{q}(n,m,k,\delta)\subseteq\mathrm{LRPC}_{q}(n,m,k,\delta^{\prime})$
if $\delta\leq\delta^{\prime}$, similarly we have that
$\text{row-}\mathrm{LRPC}_{q}(n,m,k,\rho)\subseteq\text{row-}\mathrm{LRPC}_{q}(n,m,k,\rho^{\prime})$
if $\rho\leq\rho^{\prime}$.

Let $H$ be an $(n-k)\times n$ matrix with row weight $\rho$. The weight
$\delta$ of $H$ satisfies $\rho\leq\delta\leq(n-k)\rho$. Hence, a row-LRPC
code of row weight $\rho$ may admit a parity check matrix with weight ranging
from $\rho$ up to $m$ when $(n-k)\rho\geq m$. In other words, it is clear that
$\mathrm{LRPC}_{q}(n,m,k,\rho)\subseteq\text{row-}\mathrm{LRPC}_{q}(n,m,k,\rho).$
The opposite inclusion is not true. As shown below, under certain condition it
is always possible to construct a row-LRPC code of row weight $\rho\geq 2$
that has all its parity check matrices of weight strictly greater than $\rho$.

######  Theorem 1.

Let $2\leq\rho<m$ and $2\rho-2\leq k\leq n-2$. Then we have

$\mathrm{LRPC}_{q}(n,m,k,\rho)\subsetneq\text{row-}\mathrm{LRPC}_{q}(n,m,k,\rho)$

###### Proof.

Suppose $\mathcal{A}=\langle a_{1},\ldots,a_{\rho}\rangle,\mathcal{B}=\langle
b_{1},\ldots,b_{\rho}\rangle$ are two $\rho$-dimensional
$\mathbb{F}_{q}$-vector subspace of $\mathbb{F}_{q^{m}}$ such that there
exists no $\lambda\in\mathbb{F}_{q^{m}}$ satisfying
$\mathcal{A}=\lambda\mathcal{B}$. We build an $(n-k)\times n$ parity check
matrix $H$ of row density $\rho$ in the following way

$H:=\begin{pmatrix}a_{1}&\cdots&a_{\rho}&0&\cdots&0&\mathbf{0}_{n-2\rho}\\\
0&\cdots&0&b_{1}&\cdots&b_{\rho}&\mathbf{0}_{n-2\rho}\\\
0&\cdots&0&0&\cdots&0&J\\\ \end{pmatrix},$

where $J$ is an ${(n-k-2)\times(n-2\rho)}$ matrix obtained by extending the
identity matrix of order $n-k-2$ with $k+2-2\rho$ zeros at each row, and
$\mathbf{0}_{n-2\rho}$ is a $1\times(n-2\rho)$ matrix of zeros.

Let $\mathcal{C}$ be a row-LRPC code admitting $H$ as its parity check matrix.
We can write all the possible equivalent parity check matrices for
$\mathcal{C}$ as $H^{\prime}=TH$ for some invertible $(n-k)\times(n-k)$ matrix
$T$ over $\mathbb{F}_{q^{m}}$. Since $T$ is a square invertible matrix, all
its columns have to be linearly independent. Consider the first two columns of
$T$, without loss of generality we can assume that the submatrix

$T_{\\{1,2\\},\\{1,2\\}}=\begin{pmatrix}t_{1,1}&t_{1,2}\\\ t_{2,1}&t_{2,2}\\\
\end{pmatrix}$

has non null determinant. The first $2$ rows of $TH$ can be given by

$\begin{pmatrix}t_{1,1}a_{1}&\cdots&t_{1,1}a_{\rho}&t_{1,2}b_{1}&\cdots&t_{1,2}b_{\rho}&\mathbf{x}(T)\\\
t_{2,1}a_{1}&\cdots&t_{2,1}a_{\rho}&t_{2,2}b_{1}&\cdots&t_{2,2}b_{\rho}&\mathbf{y}(T)\end{pmatrix},$

where $\mathbf{x}(T)$ and $\mathbf{y}(T)$ are two vectors of length
$(n-2\rho)$ depending on $T$ (which are both $0$-vectors in the worst case).

Let $\mathcal{T}$ be the $\mathbb{F}_{q}$-vector subspace generated by the
elements of the submatrix $TH_{\\{1,2\\},\\{1,\ldots,2\rho\\}}$ we can
describe it as

$\mathcal{T}=t_{1,1}\mathcal{A}+t_{2,1}\mathcal{A}+t_{1,2}\mathcal{B}+t_{2,2}\mathcal{B}.$

Since $T_{\\{1,2\\},\\{1,2\\}}$ is invertible, it follows that either
$t_{1,1}t_{2,2}\neq 0$ or $t_{1,2}t_{2,1}\neq 0$. We want to prove that
$\dim(\mathcal{T})>\rho$. Without loss of generality we can assume
$t_{1,1},t_{2,2}\neq 0$, since
$t_{1,1}\mathcal{A}+t_{2,2}\mathcal{B}\subseteq\mathcal{T}$ then

$\dim(t_{1,1}\mathcal{A}+t_{2,2}\mathcal{B})\leq\dim(\mathcal{T}).$

Observe that

$\displaystyle\dim(t_{1,1}\mathcal{A}+t_{2,2}\mathcal{B})$ $\displaystyle=$
$\displaystyle\dim(t_{1,1}\mathcal{A})+\dim(t_{2,2}\mathcal{B})-\dim(t_{1,1}\mathcal{A}\cap
t_{2,2}\mathcal{B})$ $\displaystyle=$ $\displaystyle
2\rho-\dim(t_{1,1}\mathcal{A}\cap t_{2,2}\mathcal{B}).$

This implies $\dim(\mathcal{T})\geq 2\rho-\dim(t_{1,1}\mathcal{A}\cap
t_{2,2}\mathcal{B})$. Since $\dim(t_{1,1}\mathcal{A}\cap
t_{2,2}\mathcal{B})\leq\rho$, then $\dim(\mathcal{T})\geq\rho$, where the
equality holds only if $\dim(t_{1,1}\mathcal{A}\cap t_{2,2}\mathcal{B})=\rho$.
This is actually equivalent to

$\mathcal{A}=t_{1,1}^{-1}t_{2,2}\mathcal{B},$

which is a contradiction since by our choice of $\mathcal{A}$ and
$\mathcal{B}$ there exists no $\lambda\in\mathbb{F}_{q^{m}}$ satisfying
$\mathcal{A}=\lambda\mathcal{B}$. Let $\mathrm{rw}(H)$ denote the rank weight
of the matrix $H$. For any invertible matrix $T$, we have

$\displaystyle\mathrm{rw}(TH)$
$\displaystyle\geq\mathrm{rw}(TH_{\\{1,2\\},\\{1,\ldots,n\\}})$
$\displaystyle\geq\mathrm{rw}(TH_{\\{1,2\\},\\{1,\ldots,2\rho\\}})>\rho.$

The desired conclusion follows. ∎

######  Remark 2.

Theorem 1 is does not cover the whole range of all parameters that make
Theorem 1’s result hold. Indeed it is quite easy to find such examples for
parameters.

Theorem 1 indicates that the extension of LRPC codes to row-LRPC codes is a
proper extension for $\rho\geq 2$. Nevertheless, this is not true for the case
when the row weight equals $1$.

######  Lemma 1.

weight1
$\small{\text{row-}\mathrm{LRPC}_{q}(n,m,k,1)=\mathrm{LRPC}_{q}(n,m,k,1).}$

###### Proof.

Let $\mathcal{C}$ be a row-LRPC code of dimension $k$ and row weight $1$. The
code $\mathcal{C}$ admits an $(n-k)\times n$ parity check matrix $H$ whose
rows have all rank weight $1$. This means $h_{i,j}\in\mathcal{H}_{i}=\langle
h_{i}\rangle$ for some $h_{i}\in\mathbb{F}_{q^{m}}^{*}$, i.e.,
$h_{i,j}=\lambda_{i,j}h_{i}$ for some $\lambda_{i,j}\in\mathbb{F}_{q}$.

The weight of the matrix $H$ is $\dim(\langle h_{1},\ldots,h_{n-k}\rangle)$
which can be any number between $1$ and $\min(n-k,m)$. We want to prove that
there exists as well another parity check matrix $H^{\prime}$ of weight $1$
for the same code $\mathcal{C}$. In fact, consider the diagonal matrix $L$
with the entries $(i,i)$ being $h_{i}^{-1}$ and all other entries zero. Then
the matrix $H^{\prime}=LH$ is still a parity check matrix of $\mathcal{C}$ and
its entries can be written as
$h_{i,j}^{\prime}=h_{i}^{-1}h_{i,j}=h_{i}^{-1}\lambda_{i,j}h_{i}=\lambda_{i,j}\in\mathbb{F}_{q}$.
Therefore the weight of $H^{\prime}$ is $1$. ∎

##  4 Decoding of Row-LRPC Codes

The syndrome decoding problem for Hamming distance is a fundamental problem in
complexity theory, which was proved to be an NP-complete problem Berlekamp,
McEliece and van Tilborg [13]. When the code is embedded with rank metric, the
corresponding rank syndrome decoding (RSD) problem is given as below.

######  Definition 5.

[6, 14] Given an $(n-k)\times n$ matrix over $\mathbb{F}_{q^{m}}$, a vector
$\mathbf{s}$ in $\mathbb{F}_{q^{m}}^{n-k}$ and an integer $r$, find a vector
$\mathbf{e}\in\mathbb{F}_{q^{m}}^{n}$ of rank weight
$\mathrm{rw}(\mathbf{e})=r$ satisfying $\mathbf{e}H^{\intercal}=\mathbf{s}.$

By embedding linear codes in the Hamming space into linear codes over an
extension field equipped with the rank metric, Gaborit and Zémor [15] proved
the hardness of the RSD problem under an unfaithful random reduction, which
indicates that if the RSD problem is in ZPP, then ZPP = NP. The RSD problem
has been intensively studied in recent years driven by its application in
post-quantum cryptosystems [16, 17, 10]. Early attacks on RSD appeared in [18,
19]. Some progress was lately made in [14], where the authors proposed two new
generic attacks: support attack and polynomial annulator attack for the RSD
problem, which outperform the combinatorial attack in [19] when the error
vector has a relatively small rank weight $r$. Very recently, algebraic
attacks with Grobner basis technique [20] and without Grobner basis technique
[21] were proposed, which have significantly advanced the RSD problem.

With the hardness of RSD problem, by far there are two general families of
rank metric codes for which (probabilistic) polynomial-time decoding
algorithms have been found. Known decoding algorithms for the Gabidulin codes
(e.g., [2]) rely on the strong algebraic structure of Gabidulin codes, and the
probabilistic decoding algorithm for the LRPC codes [22] utilizes the low-
weight property of the codes. More specifically, the decoding algorithm for
LRPC codes can be divided into two steps. The first step is to recover the
support $\mathcal{E}$ of the error vector $\mathbf{e}$ by intersecting the
vector spaces derived from the parity check matrix $H$ and the syndrome vector
$\mathbf{s}$; the second step is to reconstruct each component of the error
vector in the support $\mathcal{E}$ by re-utilizing the syndrome equation
$\mathbf{s}=\mathbf{e}H^{\intercal}$.

When a row-LRPC code of row weight $\rho\geq 2$ is considered, its parity
check matrix $H$ may have weight $\delta$ ranging from $\rho$ to
$\min(m,(n-k)\rho)$. Even in the non generic case when there exists a parity
check matrix $H^{\prime}$ for the same code of weight sufficiently small,
there is no efficient algorithm to compute $H^{\prime}$ from a randomly given
$H$. Therefore, in general, we cannot use the decoding algorithm of LRPC codes
for decoding row-LRPC codes.

In the following we shall explore the decoding process of row-LRPC codes,
particularly when the row weight is sufficiently small. Similar to the two
families of efficiently decodable rank metric codes Gabidulin codes [2] and
LRPC codes [22, 9], the critical step of decoding row-LRPC codes is to
efficiently recover the error support. Once the the error support is
recovered, it is possible to expand the $n-k$ linear equation given by the
syndrome obtaining $m(n-k)$ equation over $\mathbb{F}_{q}$ and solve the
linear system.

###  4.1 The simplest case $r=1$

The easiest case to treat is when the error support $\mathcal{E}$ has
dimension $1$. In this case, we can write
$\mathcal{E}=\langle\varepsilon\rangle$ for some
$\varepsilon\in\mathbb{F}_{q^{m}}$. The error $\mathbf{e}$ can be written as
$(\lambda_{1}\varepsilon,\ldots,\lambda_{n}\varepsilon)$ for some
$\lambda_{j}\in\mathbb{F}_{q}$. The general syndrome equation is then given by

$s_{i}=\sum_{j=1}^{n}h_{i,j}\lambda_{j}\varepsilon=\varepsilon\sum_{j=1}^{n}h_{i,j}\lambda_{j}=\varepsilon
h_{i},\qquad h_{i}\in\mathcal{H}_{i},$ (1)

where $\mathcal{H}_{i}$ is the $\rho$-dimensional $\mathbb{F}_{q}$-vector
subspace generated by the elements of the $i$-th row of $H$. From equation
(1), if $s_{i}\neq 0$ then $h_{i}\neq 0$, so we have $\varepsilon\in
s_{i}\mathcal{H}_{i}^{-1}$. That is equivalent to $\varepsilon^{-1}\in
s_{i}^{-1}\mathcal{H}$. This second formulation has an advantage that
$s_{i}^{-1}\mathcal{H}$ is an $\mathbb{F}_{q}$-vector subspace since

$\lambda_{1}s_{i}^{-1}h_{1}+\lambda_{2}s_{i}^{-1}h_{2}=s_{i}^{-1}(\lambda_{1}h_{1}+\lambda_{2}h_{2})$

belongs to $s_{i}^{-1}\mathcal{H}$ for any
$\lambda_{1},\lambda_{2}\in\mathbb{F}_{q}$, $h_{1},h_{2}\in\mathcal{H}$. By
intersecting several sets of the form $s_{i}^{-1}\mathcal{H}_{i}$, we will
likely obtain the $\mathbb{F}_{q}$-vector space
$\mathbb{F}_{q}\varepsilon^{-1}$ from which it’s easy to find the support
$\mathcal{E}$.

###  4.2 Cramer’s rule decoding for $r=2$

We will present a method to recover the error support when the error has rank
$r\geq 2$. Before treating the general case, we shall start with the case
$r=2$ as a simpler illustration in this subsection.

Notice that, as previously shown, the set
$s\mathcal{H}=\\{sh|h\in\mathcal{H}\\}$ is a $\mathbb{F}_{q}$-vector subspace
if $\mathcal{H}$ is a $\mathbb{F}_{q}$-vector subspace. Moreover, it is clear
that $s\mathcal{H}=\langle sh_{1},\ldots,sh_{\rho}\rangle$ when
$\mathcal{H}=\langle h_{1},\ldots,h_{\rho}\rangle$ since any element in
$s\mathcal{H}$ can be linearly expressed by $sh_{1},\ldots,sh_{\rho}$. This
fact will be frequently used in the subsequent decoding process.

Let $\mathcal{E}=\langle\varepsilon_{1},\varepsilon_{2}\rangle$ be the error
support of $\mathbf{e}$ and let $H$ be a parity check matrix of row weight
$\rho$. The general syndrome equation will be of the form

$h_{i,1}\varepsilon_{1}+h_{i,2}\varepsilon_{2}=s_{i},\qquad
h_{i,1},h_{i,2}\in\mathcal{H}_{i}$

where $\mathcal{H}_{i}$ is the $\rho$-dimensional $\mathbb{F}_{q}$-vector
subspace generated by the elements of the $i$-th row of the parity check
matrix $H$. When $s_{i}\neq 0$ we can multiply the whole equation by
$s_{i}^{-1}$. Let $\mathcal{A}_{i}$ be an $\mathbb{F}_{q}$-vector space given
by

$\mathcal{A}_{i}=\begin{cases}s_{i}^{-1}\mathcal{H}_{i},&\text{if }s_{i}\neq
0,\\\ \mathcal{H}_{i},&\text{otherwise}.\end{cases}$

In this way, we can rewrite the equation
$H\mathbf{e}^{\intercal}=\mathbf{s}^{\intercal}$ as follows

$A\cdot\mathbf{e}^{\intercal}=\mathbf{v}\Leftrightarrow\begin{pmatrix}a_{1,1}&a_{1,2}\\\
a_{2,1}&a_{2,2}\\\ \vdots&\vdots\\\
a_{n-k,1}&a_{n-k,r}\end{pmatrix}\begin{pmatrix}\varepsilon_{1}\\\
\varepsilon_{2}\end{pmatrix}=\begin{pmatrix}v_{1}\\\ v_{2}\\\ \vdots\\\
v_{n-k}\end{pmatrix},$ (2)

where $a_{i,j}\in\mathcal{A}_{i}$ and
$\mathbf{v}=(v_{1},\ldots,v_{n-k})^{\intercal}$ is a vector over
$\mathbb{F}_{q^{m}}$ such that $v_{i}=0$ when $s_{i}=0$ and is $1$ otherwise.

Suppose that $v_{i}=v_{j}=1$, from (2) we can extract the linear system

$A_{\\{i,j\\},\\{1,2\\}}\cdot\mathbf{e}^{\intercal}=\begin{pmatrix}a_{i,1}&a_{i,2}\\\
a_{j,1}&a_{j,2}\\\ \end{pmatrix}\begin{pmatrix}\varepsilon_{1}\\\
\varepsilon_{2}\end{pmatrix}=\begin{pmatrix}1\\\ 1\end{pmatrix}.$ (3)

There are at most $q^{4\rho}$ possible configuration of the matrix
$A_{\\{i,j\\},\\{1,2\\}}$. If the matrix $A_{\\{i,j\\},\\{1,2\\}}$ is
invertible, we can solve (3) by Cramer’s rule [23] as

$\varepsilon_{1}=\frac{\det\begin{pmatrix}1&a_{i,2}\\\
1&a_{j,2}\end{pmatrix}}{\det\begin{pmatrix}a_{i,1}&a_{i,2}\\\
a_{j,1}&a_{j,2}\end{pmatrix}},\quad\varepsilon_{2}=\frac{\det\begin{pmatrix}a_{i,1}&1\\\
a_{j,1}&1\end{pmatrix}}{\det\begin{pmatrix}a_{i,1}&a_{i,2}\\\
a_{j,1}&a_{j,2}\end{pmatrix}}.$ (4)

After calculation of determinants this reduces to

$\varepsilon_{1}=\frac{a_{i,2}-a_{j,2}}{a_{i,1}a_{j,2}-a_{i,2}a_{j,1}},\quad\varepsilon_{2}=\frac{a_{i,1}-a_{j,1}}{a_{i,1}a_{j,2}-a_{i,2}a_{j,1}}.$
(5)

We define the set

$B_{i,j}=\left\\{\frac{a_{i,2}-a_{j,2}}{a_{i,1}a_{j,2}-a_{i,2}a_{j,1}}\,:\,a_{i,*}\in\mathcal{A}_{i},\,a_{j,*}\in\mathcal{A}_{j}\right\\}.$
(6)

Recall that $\mathcal{A}.\mathcal{B}$ denotes the product vector space
generated by the products of elements in the vector spaces $\mathcal{A}$ and
$\mathcal{B}$. For the set $B_{i,j}$ in (6), the denominators in the division
satisfy

$a_{i,1}a_{j,2}-a_{i,2}a_{j,1}\in\mathcal{A}_{i}.\mathcal{A}_{j}$

since $a_{i,*}a_{j,*}\in\mathcal{A}_{i}.\mathcal{A}_{j}$. This leads to an
observation that

$B_{i,j}\subseteq\frac{\mathcal{A}_{i}+\mathcal{A}_{j}}{\mathcal{A}_{i}.\mathcal{A}_{j}}.$

In addition, due to the symmetry of the equations in (5), we have

$\mathcal{E}\subsetneq
B_{i,j}\subseteq\frac{\mathcal{A}_{i}+\mathcal{A}_{j}}{\mathcal{A}_{i}.\mathcal{A}_{j}}.$

Since this is true for any $B_{i,j}$ we can intersect different $B_{i,j}$ in
order to get a smaller set where to look for the error support $\mathcal{E}$.
Since the sets $B_{i,j}$ depend on $4$ parameters in $\mathbb{F}_{q}$-vector
spaces of dimension $\rho$, their size is upper bounded by $q^{4\rho}$.

Following the method just explained, we run some MAGMA simulations for the
case $r=2$, $\rho=2$, $n=m=20$ and $q=2,3$. Using $4$ rows in the parity check
matrix it was possible to recover the error support in $707/1000$ cases for
$q=2$ and in $954/1000$ cases for $q=3$. Increasing the number of parity check
rows to $6$, the results improved to $988/1000$ and $998/1000$, respectively.
In all of the above mentioned cases, it was possible to recover the error
support after the intersection of just two $B_{i,j}$ sets.

In Section 5 we will explain why and how increasing the number of rows of the
parity check matrix can help us achieve higher success rates in decoding.

###  4.3 Cramer’s rule for the general case

The decoding process with Cramer’s rule for $r=2$ can be easily extended to
general cases.

Let $\mathcal{E}=\langle\varepsilon_{1},\ldots,\varepsilon_{r}\rangle$ be the
support of the error $\mathbf{e}=(e_{1},\ldots,e_{n})$. Given $H$ a parity
check matrix of row weight $\rho$ we can rewrite the $i$-th syndrome equation
$\mathbf{e}\cdot\mathbf{h}_{i}^{\intercal}=s_{i}$ as

$\varepsilon_{1}h_{i,1}+\cdots+\varepsilon_{r}h_{i,r}=s_{i}\qquad
h_{i,*}\in\mathcal{H}_{i},$ (7)

where $\mathcal{H}_{i}=\langle\mathbf{h}_{i}\rangle$ is the $\rho$-dimensional
vector space generated by the $i$-th row of $H$.

As we did in the previous section, we can manipulate each equations
multiplying on both sides by $s_{i}^{-1}$ when $s_{i}\neq 0$. We rename the
variables and the $\mathbb{F}_{q}$-vector subspace with the same notation. In
this way we get a system that can be described by the matrix equation
$A\cdot\mathbf{e}^{\intercal}=\mathbf{v}$, more specifically,

$\begin{pmatrix}a_{1,1}&\cdots&a_{1,r}\\\ a_{2,1}&\cdots&a_{2,r}\\\
\vdots&\ddots&\vdots\\\
a_{n-k,1}&\cdots&a_{n-k,r}\end{pmatrix}\begin{pmatrix}\varepsilon_{1}\\\
\vdots\\\ \varepsilon_{r}\end{pmatrix}=\begin{pmatrix}v_{1}\\\ v_{2}\\\
\vdots\\\ v_{n-k}\end{pmatrix}.$ (8)

This time we can select $r$ rows of (8) and solve by the Cramer’s rule. Let
$I=\\{i_{1},\ldots,i_{r}\\}\subseteq\\{1,\ldots,n-k\\}$ a set of $r$ indices.
The Cramer’s rule for the sub-system given by the the rows indexed by $I$ is

$\varepsilon_{j}=\frac{\det(A_{I,\\{1,\ldots,j-1\\}}|\mathbf{v}|A_{I,\\{j+1,\ldots,r\\}})}{\det(A_{I,\\{1,\ldots,r\\}})},$
(9)

where the matrix $A_{I,\\{1,\ldots
j-1\\}}|\mathbf{v}|A_{I,\\{j+1,\ldots,r\\}}$ is the matrix
$A_{I,\\{1,\ldots,r\\}}$ where the $j$-th column is substituted by
$\mathbf{v}$. In this way we can define the set
$B_{I}\subseteq\mathbb{F}_{q^{m}}$ of size upper bounded by $q^{r^{2}\rho}$
as:

$B_{I}=\left\\{\frac{\det(A_{I,\\{1,\ldots,r-1\\}}|\mathbf{v})}{\det(A_{I,\\{1,\ldots,r\\}})}\bigg{|}\quad
a_{i,j}\in\mathcal{A}_{i}\right\\}$ (10)

where $a_{i,j}$ indicates the entries of $A_{I,\\{1,\ldots,r\\}}$ at the row
$i$ and column $j$. Notice that since the variables $a_{i,j}$ are uniformly
chosen from all the possible values of the $\mathbb{F}_{q}$-vector space
$\mathcal{A}_{i}$, the quantity $\det(A^{{}^{\prime}}_{I,\\{1,\ldots
j-1\\}}|\mathbf{v}|A^{{}^{\prime}}_{I,\\{j+1,\ldots,r\\}})$ with
$j=1,\ldots,r$ can be obtained as well as
$\det(A_{I,\\{1,\ldots,r-1\\}}|\mathbf{v})$, where $A$ is obtained from
$A^{\prime}$ by renaming the columns after $j$ and multiplying one column by
$-1$ when $\mathbb{F}_{q}$ has odd characteristic.

Each $B_{I}$ contains the error support $\mathcal{E}$. We can then intersect
some different $B_{I}$’s to obtain a smaller set that can be exactly
$\mathcal{E}$ or at least a smaller set where we know $\mathcal{E}$ is
contained. The situation where we expect this algorithm to work is when the
cardinality of $B_{I}$ is much smaller than $\mathbb{F}_{q^{m}}$ (i.e.
$r^{2}\rho\ll m$).

If we intersect two sets $B_{I}$ and $B_{J}$ such that $I\cap J\neq\emptyset$,
we can expect to have a bigger intersection. Suppose without loss of
generality that the index $1$ belongs to both $I$ and $J$. Among all the
possible configurations of $A_{I,\\{1,\ldots,r\\}}$ and
$A_{J,\\{1,\ldots,r\\}}$ we will have one where the first row of both matrices
is $(0,\ldots,0,a_{1,r})$. When we apply the Cramer’s rule to that matrix we
will obtain $a_{1,r}^{-1}$, then we know without the need of computing that
$\mathcal{A}_{1}^{-1}\subset B_{I}\cap B_{J}$.

For this reason, in order to do as few intersections as possible, it is better
to use disjointed indexes. If we want to have at least two disjointed sets of
indexes, we need $n-k\geq 2r$. Suppose we managed to recover the error
support, to recover each component of the error $\mathbf{e}$ we need to solve
an expanded linear system of $nr$ unknowns and $m(n-k)$ equations. Therefore
we need that $(n-k)\geq\lceil\frac{nr}{m}\rceil$. Then we have the condition
on the dimension $k$

$k\leq\min\left(n-2r,n-\left\lceil\frac{nr}{m}\right\rceil\right).$

The hardest part of this algorithm is ruled by the error support recovery. In
particular it is given by the cost of computing the sets $B_{I}$, to obtain
each element of $B_{I}$ we need to solve a linear system in $r$ unknowns. We
can then give an estimation of the total complexity in the order of
$O(r^{3}q^{r^{2}\rho})$.

##  5 Decoding failures

We implement the algorithm described in the previous section on MAGMA on an
easy instance of the problem. In particular we created different instances of
row-LRPC codes of row weight $2$ and we did correct errors of rank weight $2$.
During the error support recovery, we noticed two possible causes of failure.

###  5.1 Syndrome contains some zeroes

Let $\mathcal{E}=\langle\varepsilon_{1},\ldots,\varepsilon_{r}\rangle$ be the
error support. The generic equation of the system is

$a_{i,1}\varepsilon_{1}+\cdots+a_{i,r}\varepsilon_{r}=s_{i}.$ (11)

We can consider $a_{i,*}$ as independent variables uniformly distributed over
$\mathcal{H}_{i}$, then $a_{i,j}\varepsilon_{j}$ are independent uniformly
distributed over $\mathcal{H}_{i}\varepsilon_{j}$. Therefore $s_{i}$ is a
variable uniformly distributed over the product space
$\mathcal{H}_{i}.\mathcal{E}$ which, with an high probability, has dimension
$r\rho$ when $r\rho\ll m$ [24].
By our hypothesis
$\mathcal{H}_{i}=\langle\alpha_{1},\ldots,\alpha_{\rho}\rangle$ then each
$a_{i,j}$ can be rewritten as
$\sum_{i\in\\{1,\ldots,\rho\\}}\lambda_{i,j}\alpha_{i}$ for some
$\lambda_{i,j}\in\mathbb{F}_{q}$. The equation (11) can be rewritten as

$\sum_{j\in\\{1,\ldots,r\\}}\sum_{l\in\\{1,\ldots,\rho\\}}\lambda_{l,j}\alpha_{l}\varepsilon_{j}=s_{i}\quad\lambda_{l,j}\in\mathbb{F}_{q}.$

Suppose now that $s_{i}=0$ and $r\rho\ll m$, with an high probability,
$\alpha_{l}\varepsilon_{j}$ are all $\mathbb{F}_{q}$-linearly independent
which implies $\lambda_{i,j}=0\quad\forall
i\in\\{1,\ldots,\rho\\},\quad\forall j\in\\{1,\ldots,r\\}$, therefore
$a_{i,1}=\cdots=a_{i,\rho}=0$.
Unfortunately in the case $s_{i}=0$ we cannot get any information about the
error support. Indeed, if all the coefficients $a_{i,*}$ are zero, the
equation (11) is trivially satisfied by all the possible sets of $r$ elements
in $\mathbb{F}_{q^{m}}$.

If we consider $s_{i}$ as independent uniformly distributed variables over the
$\mathbb{F}_{q}$-vector subspace $\mathcal{H}_{i}.\mathcal{E}$ of dimension
$r\rho$, we have the probability $P(s_{i}=0)=q^{-r\rho}$. We expect to have at
least $t$ non-zero components of the syndrome with probability
$\sum_{i=t}^{n-k}\binom{n-k}{i}(1-q^{-r\rho})^{i}q^{-r\rho(n-k-i)}$.

The equations leading to zero will not be totally discarded as they can still
be used, after the recovery of the error support, to help us recovering the
actual error component by component.

###  5.2 Intersecting subspaces

For simplicity here we will analyze the case when the error support has rank
$2$. In the error support recovery, after we discard all the equation
corresponding to $s_{i}=0$, there is still a possible issue we will show in
this paragraph.

Let $\mathcal{E}=\langle\varepsilon_{1},\varepsilon_{2}\rangle$ be the error
support, consider the generic system:

$\begin{cases}a_{i,1}\varepsilon_{1}+a_{i,2}\varepsilon_{r}=1\\\
a_{j,1}\varepsilon_{1}+a_{j,2}\varepsilon_{r}=1\end{cases}$ (12)

where $a_{i,*}\in\mathcal{A}_{i}=\mathcal{H}_{i}s_{i}^{-1}$ and
$a_{j,*}\in\mathcal{A}_{j}=\mathcal{H}_{j}s_{j}^{-1}$.

Suppose that $\mathcal{A}_{i}$ and $\mathcal{A}_{j}$ intersect non-trivially,
the associated set $B_{i,j}$ will be smaller than what we usually expect. A
smaller set is not necessarily a good news. At the contrary, what we observed
with MAGMA, is that many times this set does contain only part of the error
support. When we intersect such a set with the other sets containing the whole
support, we recover only part of the support instead of the whole as we wish.

Given two randomly generated $\mathbb{F}_{q}$-vector subspaces of dimension
$\rho\ll m$, the probability of having a non-trivial intersection between the
two should be negligible. What we observe with MAGMA, suggests that in our
case this happens with a much higher frequency than one could naively expect.

Consider better the system (12), without loss of generality assume $i=1$ and
$j=2$

$\begin{cases}a_{1,1}\varepsilon_{1}+a_{1,2}\varepsilon_{2}=1\\\
a_{2,1}\varepsilon_{1}+a_{2,2}\varepsilon_{2}=1\end{cases}$

Consider the case $a_{1,2}=a_{2,2}=0$, this would imply
$a_{1,1}\varepsilon_{1}=a_{2,1}\varepsilon_{1}=1$, then
$\varepsilon_{1}^{-1}\in\mathcal{A}_{1}\cap\mathcal{A}_{2}$. A similar
conclusion can be made from $a_{1,1}=a_{2,1}=0$. A bit less straight forward
is the case $a_{1,1}=a_{1,2}\neq 0$. In this case we have
$(\varepsilon_{1}+\varepsilon_{2})^{-1}\in\mathcal{A}_{1}$, if we also have
$a_{2,1}=a_{2,2}\neq 0$ then
$(\varepsilon_{1}+\varepsilon_{2})^{-1}\in\mathcal{A}_{1}\cap\mathcal{A}_{2}$.
The coefficients $a_{i,*}$ all belong to $\mathcal{A}_{i}$ which, by our
hypothesis, is a small set of size $q^{\rho}$. This explains why we observe
non trivial intersections so often.

In the general case, after we consider a system of $r$ rows, we will obtain an
$r\times r$ matrix $A$ whose entries belong to some $\mathbb{F}_{q}$-vector
subspace of dimension $\rho$, that multiplies the vector
$\mathbf{e}=(\varepsilon_{1},\ldots,\varepsilon_{r})$ containing a base of the
error support. Intuitively, in the case the $j$-th column of $A$ is all $0$s,
we have no hope to get some information about $\varepsilon_{j}$. The
likelihood of this event decreases exponentially in $r$. In the case $r>2$, we
cannot find a relation between the $\mathbb{F}_{q}$-vector subspaces
$\mathcal{A}_{i}$ that tells us something about the matrix $A$. Increasing the
number of equation will increase the probability of success. In the case of
intersecting spaces, it will give us more $\mathbb{F}_{q}$-vector subspaces
among which we can find more disjointed couples. In the case of too many
zeroes in the syndrome, increasing the number of equations will increase the
probability of having enough non-zero components.

##  6 Conclusion

This paper introduces a new family of rank metric codes, which can be seen as
an extended variant of LRPC codes introduced in [22]. We showed that the
proposed row-LRPC codes contain LRPC codes as a strictly proper subfamily when
row weight is no less than $2$ under certain condition, which confirms the
validity of the extension. The decoding of row-LRPC codes is also studied by
intersecting certain sets obtained from the syndrome equation, which has
complexity upper bounded by $O(r^{3}q^{r^{2}\rho})$, where $\rho$ is the row
weight of the code in question and $r$ is the rank weight of the error vector.

Due to high complexity of the current decoding approach, we will continue
exploring the possibility of a polynomial-time decoding algorithm for the row-
LRPC codes, which by then would have interesting applications in
cryptography..

## References

  * [1] P. Delsarte, “Bilinear forms over a finite field, with applications to coding theory,” _Journal of Combinatorial Theory, Series A_ , vol. 25, no. 3, pp. 226 – 241, 1978.
  * [2] E. M. Gabidulin, “Theory of codes with maximum rank distance,” _Problemy Peredachi Informatsii_ , vol. 21, no. 1, pp. 3–16, 1985.
  * [3] R. M. Roth, “Maximum-rank array codes and their application to crisscross error correction,” _IEEE Transactions on Information Theory_ , vol. 37, no. 2, pp. 328–336, 1991.
  * [4] D. Silva, F. R. Kschischang, and R. Koetter, “A rank-metric approach to error control in random network coding,” _IEEE Transactions on Information Theory_ , vol. 54, no. 9, pp. 3951–3967, Sept 2008.
  * [5] P. Lusina, E. Gabidulin, and M. Bossert, “Maximum rank distance codes as space-time codes,” _IEEE Transactions on Information Theory_ , vol. 49, no. 10, pp. 2757–2760, 2003.
  * [6] E. M. Gabidulin, A. V. Paramonov, and O. V. Tretjakov, “Ideals over a non-commutative ring and their application in cryptology,” in _Advances in Cryptology – EUROCRYPT’91_ , D. W. Davies, Ed. Springer, 1991, pp. 482–489.
  * [7] R. Gallager, “Low-density parity-check codes,” _IRE Transactions on Information Theory_ , vol. 8, no. 1, pp. 21–28, 1962.
  * [8] P. Gaborit, A. Hauteville, D. H. Phan, and J.-P. Tillich, “Identity-based encryption from codes with rank metric,” in _Advances in Cryptology – CRYPTO 2017_ , J. Katz and H. Shacham, Eds. Springer International Publishing, 2017, pp. 194–224.
  * [9] N. Aragon, P. Gaborit, A. Hauteville, O. Ruatta, and G. Zémor, “Low rank parity check codes: New decoding algorithms and applications to cryptography,” _IEEE Transactions on Information Theory_ , vol. 65, no. 12, pp. 7697–7717, 2019.
  * [10] C. A. Melchor, N. Aragon, M. Bardet, S. Bettaieb, L. Bidoux, O. Blazy, J.-C. Deneuville, P. Gaborit, A. Hauteville, A. Otmani, O. Ruatta, J.-P. Tillich, and G. Zémor, “ROLLO (merger of Rank-Ouroboros, LAKE and LOCKER),” in _Second round submission to the NIST post-quantum cryptography call_ , April, 2020.
  * [11] P. Gaborit, O. Ruatta, J. Schrek, and G. Zémor, “Ranksign: an efficient signature algorithm based on the rank metric,” in _Post-Quantum Cryptography_ , M. Mosca, Ed. Springer International Publishing, 2014, pp. 88–107.
  * [12] B. Biswas and N. Sendrier, “McEliece cryptosystem implementation: theory and practice,” in _Post-Quantum Cryptography_ , J. Buchmann and J. Ding, Eds. Berlin, Heidelberg: Springer Berlin Heidelberg, 2008, pp. 47–62.
  * [13] E. Berlekamp, R. McEliece, and H. van Tilborg, “On the inherent intractability of certain coding problems (corresp.),” _IEEE Transactions on Information Theory_ , vol. 24, no. 3, pp. 384–386, 1978.
  * [14] P. Gaborit, O. Ruatta, and J. Schrek, “On the complexity of the rank syndrome decoding problem,” _IEEE Transactions on Information Theory_ , vol. 62, no. 2, pp. 1006–1019, 2016.
  * [15] P. Gaborit and G. Zémor, “On the hardness of the decoding and the minimum distance problems for rank codes,” _IEEE Transactions on Information Theory_ , vol. 62, no. 12, pp. 7245–7252, 2016.
  * [16] P. Gaborit, O. Ruatta, J. Schrek, and G. Zémor, “New results for rank-based cryptography,” in _Progress in Cryptology – AFRICACRYPT 2014_ , D. Pointcheval and D. Vergnaud, Eds. Springer International Publishing, 2014, pp. 1–12.
  * [17] P. Loidreau, “A new rank metric codes based encryption scheme,” in _Post-Quantum Cryptography_ , T. Lange and T. Takagi, Eds. Springer International Publishing, 2017, pp. 3–17.
  * [18] F. Chabaud and J. Stern, “The cryptographic security of the syndrome decoding problem for rank distance codes,” in _Advances in Cryptology — ASIACRYPT ’96_ , K. Kim and T. Matsumoto, Eds. Berlin, Heidelberg: Springer Berlin Heidelberg, 1996, pp. 368–381.
  * [19] A. V. Ourivski and T. Johansson, “New technique for decoding codes in the rank metric and its cryptography applications,” _Probl. Inf. Transm._ , vol. 38, no. 3, pp. 237–246, 2002.
  * [20] M. Bardet, P. Briaud, M. Bros, P. Gaborit, V. Neiger, O. Ruatta, and J. Tillich, “An algebraic attack on rank metric code-based cryptosystems,” in _Advances in Cryptology - EUROCRYPT 2020_ , ser. Lecture Notes in Computer Science, A. Canteaut and Y. Ishai, Eds., vol. 12107. Springer, 2020, pp. 64–93.
  * [21] M. Bardet, M. Bros, D. Cabarcas, P. Gaborit, R. A. Perlner, D. Smith-Tone, J. Tillich, and J. A. Verbel, “Improvements of algebraic attacks for solving the rank decoding and minrank problems,” in _Advances in Cryptology - ASIACRYPT 2020_ , ser. Lecture Notes in Computer Science, S. Moriai and H. Wang, Eds., vol. 12491. Springer, 2020, pp. 507–536.
  * [22] Philippe Gaborit, Gaétan Murat, Olivier Ruatta, and Gilles Zémor, “Low rank parity check codes and their application to cryptography. in proceedings of the workshop on coding and cryptography WCC’2013 Bergen Norway 2013. available on www.selmer.uib.no/wcc2013/pdfs/gaborit.pdf.”
  * [23] A. A. Kosinski, “Cramer’s rule is due to Cramer,” _Mathematics Magazine_ , vol. 74, no. 4, pp. 310–312, 2001.
  * [24] N. Aragon, P. Gaborit, A. Hauteville, and J. Tillich, “A new algorithm for solving the rank syndrome decoding problem,” in _2018 IEEE International Symposium on Information Theory (ISIT)_ , 2018, pp. 2421–2425.


# Detecting directional forces in the evolution of grammar: A case study of
the English perfect with intransitives across EEBO, COHA, and Google Books

Shimpei Okuda  Graduate School of Informatics, Nagoya University, Furo-cho,
Chikusa-ku, Nagoya 464-8601, Japan  Michio Hosaka∗ Department of English
Language and Literature, Nihon University, 3-25-40 Sakurajyosui, Setagaya-ku,
Tokyo, Japan  Kazutoshi Sasahara† School of Environment and Society, Tokyo
Institute of Technology, 3-3-6 Shibaura Minato-ku, Tokyo 108-0023, Japan

(Corresponding authors:
∗<EMAIL_ADDRESS>
†<EMAIL_ADDRESS>
###### Abstract

Languages have diverse characteristics that have emerged through evolution. In
modern English grammar, the perfect is formed with have+PP (past participle),
but in earlier English the be+PP form also existed. It is widely recognised
that the auxiliary verb BE was replaced by HAVE throughout evolution, except
for some special cases. However, whether this evolution was caused by natural
selection or random drift is still unclear. Here we examined directional
forces in the evolution of the English perfect with intransitives by combining
three large-scale data sources: EEBO (Early English Books Online), COHA
(Corpus of Historical American English), and Google Books. We found that most
intransitive verbs exhibited an apparent transition from be+PP to have+PP,
most of which were classified as ‘selection’ by a deep neural network-based
model. These results suggest that the English perfect could have evolved
through natural selection rather than random drift, and provide insights into
the cultural evolution of grammar.

### Keywords

cultural evolution, evolutionary force, human language, natural selection,
random drift

##  1 Introduction

###  1.1 Cultural evolution of human language

There are approximately 7000 languages in the world today, and these languages
are unique in terms of their phonology, morphology, and grammar. Such
diversity is typically observed in the biological world. Organisms began as
prokaryotes, and through the processes of ‘macro-evolution’ (major changes in
species occurring over a long period) and ‘micro-evolution’ (small changes in
diversity within species), these processes resulted in species of various
creatures. In evolutionary linguistics, it is believed that modern languages
have been diversified through a similar process [1]; that is, language ability
could have emerged in the human species through the macro-evolution, and then
it could have led to diverse modern languages in the world through the micro-
evolution. The former is called ‘biological evolution’, and the latter is
called ‘cultural evolution’. This paper focuses on the cultural evolution of
language, especially grammar.

The study of cultural evolution has dealt with various subjects, including
institutions, morality, and religion [2]. Unlike the evolution of biological
species, most of these phenomena are not visible in traits. It is often said
that ‘language does not fossilise’ [3]. The cultural evolution of language has
been particularly difficult to study. However, the progress of computational
models, analytical methods, and laboratory experiments has led to the
development of research on the cultural evolution of language [4, 5, 6]. In
addition, the appearance of big data and information technology has allowed us
to quantify the dynamics of the cultural evolution of language [7], as already
seen in the evolution of biological species. Notably, Newberry et al. [8]
analysed the historical changes in the inflexions of verbs and the auxiliary
verb DO in English, demonstrating that these changes are caused not only by
natural selection but also by random drift (neutral evolution).

In this paper, we focus on a specific case of the cultural evolution of
grammar—auxiliary verb selection in the evolution of the English perfect.
Although this phenomenon was previously documented with the analysis of small
corpora, we quantify it by using three large-scale data sources, providing
insights into the evolutionary dynamics of English grammar.

###  1.2 Auxiliary verb selection and evolution of the perfect

It is widely recognised that, in most cases, the auxiliary verb BE was
replaced by HAVE throughout the evolution of the English perfect. This
phenomenon is referred to as ‘auxiliary selection’.

First, we review previous studies of the perfect in historical linguistics. In
many languages, notably Indo-European languages, the auxiliary verb BE or HAVE
is used to construct the perfect [9, 10, 11]. In such cases, certain
restrictions can be observed in which an auxiliary verb is chosen. In general,
HAVE is chosen for transitive verbs, while BE is chosen for intransitive verbs
(see (1)).

(1) | a. |  Ria heeft de schuur geverfd. | (Dutch)
---|---|---|---
|  | Ria has the shed painted. |
|  | ‘Ria has painted the shed.’ |
| b. |  Onze nieuwe piano is eindelijk gearriveerd. | (Dutch)
|  | our new piano is finally arrived |
|  | ‘Our new piano has finally arrived’ | (Ackema and Sorace 2017:2)

However, in the case of intransitive verbs, BE is selected for the
unaccusative verb, which assigns a theme role to an underlying object (see
(2)), while HAVE is selected for the unergative verb which assigns an agent
role to its subject (see (3)). This phenomenon is referred to as auxiliary
selection.

(2) | a. | Ma sœur est arrivée/*a arrivé en retard. | (French)
---|---|---|---
|  | my sister is arrived/has arrived late |
| b. | Der Zug ist/*hat spät angekommen. | (German)
|  | the train is/has late arrived | (Sorace 2004: 256)
(3) | a. | Les ouvriers ont travaillé/*sont travaillés toute la nuit. | (French)
|  | the workmen have worked/are worked whole the night |
| b. | Kurt hat/*ist den ganzen Tag gearbeitet. | (German)
|  | Kurt has/is the whole day worked |  (ditto)

In English, however, there has been a substantial change in the auxiliary verb
selection during the process of cultural evolution [12, 13]. In modern
English, only HAVE is used as an auxiliary verb (see (4)), but the perfect use
of the auxiliary verb BE was often seen (see (5)) from Old English through
Late Modern English. It is recognised today that the number of such perfect
variations using the verb BE decreased markedly in the 19th century [14].

(4) | a. |  John has/*is eaten pizza. |
---|---|---|---
| b. |  John has/*is worked for an hour. |
| c. |  John has/*is arrived | (McFadden 2007: 675)
(5) | a. |  oþþæt wintra bið þusend urnen |
|  | until winters(GEN) is thousand run |
|  | ‘until a thousand years have passed’ | (Phoen 363; Denison 1993: 359)
| b. |  Whanne he escaped was |  (Chaucer, CT.Mk. VII.2735; ditto)
| c. |  yet Benedicke was such another; and now is he become |
|  | a man. |  (Shakespeare, Ado III.iv.86; ditto)

To analyse the frequency and structural changes of the HAVE and BE perfect,
McFadden and Alexiadou [15] used corpora of Old English through Early Modern
English (The York-Toronto-Helsinki Parsed Corpus of Old English Prose (YCOE),
The Penn-Helsinki Parsed Corpus of Middle English (PPCME2), The Penn-Helsinki
Parsed Corpus of Early Modern English (PPCEME), and McFadden [16] used a
corpus spanning Early Modern to Modern English (PPCMBE). These corpora are
relatively small: YCOE (1.5 M), PPCME2 (1.2 M), PPCEME (1.7 M), and PPCMBE
(1.4 M). As a result of these analyses, it was observed that the total
frequency of the HAVE perfect increased in Late Middle English (1420–1500),
and that of the BE perfect did not decrease until around 1700 and then
decreased in the 19th century (1810–1861).111McFadden and Alexiadou (2010:
422) explains the time gap between the rise of HAVE-perfect and the fall of
BE-perfect as follows: “1350 was not actually the start of the loss of be, but
the start of an expansion of have at the expense of the simple past. Indeed,
we have given quantitative evidence that the frequency of the be periphrasis
was stable throughout Middle English and Early Modern English, that is, up to
around 1700.”

Previous studies succeeded in identifying some of the factors that influenced
the increase of the HAVE perfect, which include linguistic factors such as
past counterfactual modality, iterative, durative and atelic meanings, as well
as telic eventualities and extralinguistic factors such as chronology and text
type [12, 14, 15, 16]. Here, we add quantitative evidence of directional
forces for grammar evolution, providing detailed dynamics of each single verb
(i.e., evolutionary speed and trajectory) based on computational analysis of
large-scale data sources.

With reference to these studies, we aimed to assess whether this change was
caused by natural selection or random drift using a large-scale dataset. We
combined three large-scale data sources (as explained in the next section) to
detect evolutionary forces from systematically selected verbs in the English
perfect using a robust method. The results complement previous findings and
provide insights into the cultural evolution of English grammar.

##  2 Data and methods

###  2.1 Data

Table 1 summarises the large-scale data sources used for this study. The first
is Early English Books Online (EEBO), a historical corpus that contains
British-related printed matters published between the 15th and 17th centuries
[17]. EEBO contains various texts spanning genres, such as literature,
history, philosophy, science, politics, law, and economy. Given a query, the
EEBO database returns matched sentences with the relevant information,
including the corresponding document IDs, pages, and years, from which we can
compute the frequency of matched cases (counts).

The second is the Corpus of Historical American English (COHA), comprising
American English texts in fiction, popular magazines, newspapers, and non-
fiction (books) spanning the 19th to 21st centuries [18].222COHA was updated
in 2021. We used the 2010 version. COHA is one of the largest structured
historical corpora of American English, and it enables performing queries and
providing the matched sentences with their time stamps. Similarly, from this,
we can compute the frequency of matched cases (counts).

The third is $N$-gram data from Google Books Ngram Viewer site (English, 2012
version) [19] that is used to incorporate data spanning the 18th to 21st
centuries, which EEBO and COHA do not cover. In this paper, we simply refer to
this as Google Books data. Using web scraping, we searched and collected both
British and American English publications between 1700 and 2000 at the site,
as EEBO covers British English 1473–1700 and COHA does American English
1810–2009.

It should be noted that the collected Google Books Ngram Viewer site does not
have the Part-of-Speech (POS) tag function; given a simple query, it only
provides the relative frequency of words/phrases (NOT raw counts), without
matched source texts. Thus, we need to adjust the resulting frequencies of
Google Books $N$-gram search based on those from the two other corpora, which
we call ‘scaling’ and explain later.

Table 1: Three data sources used in this study. Corpus | Year | Size (words)
---|---|---
1\. Early English Books Online (EEBO) | 1473–1700 | 7.55 M
2\. Corpus of Historical American English (COHA) | 1810–2009 | 400 M
3\. Google Books | 1700–2000 | 468 B

###  2.2 Intransitive verbs

As the be+PP (past participle) construction is used not only in the perfect
but also in the passive, it is difficult for computer programs to
automatically distinguish between these in large-scale corpora (EEBO and COHA)
and especially in big data without grammatical tags (Google Books). Therefore,
in this study, we focused only on ‘intransitive verbs’ as targets to
accurately detect the perfect construction of be+PP. This is because
intransitive verbs are not used in the passive, in principle. Due to this
restriction, we cannot cover diverse expressions of the perfect, but instead
we can accurately count the occurrences of the perfect in the datasets. To
select intransitive verbs, we used the Longman Dictionary of Contemporary
English (LDOCE) Online [20], which is a standard dictionary of the English
language.

###  2.3 Selection of target verbs

To quantify the cultural evolution of the English perfect, we systematically
selected intransitive verbs that appeared in all three datasets and had a high
frequency within them. In this study, we used two groups of target verbs: (A)
verbs selected based on the frequency and (B) verbs that did not meet the
conditions in (A) but were used in prior studies, as shown next.

#### (A) Verb group selected based on frequency

The Corpus of Contemporary American English (COCA; 400M tokens between 1990
and 2019) provides a full list of the most frequent contemporary words
($n=60,000$), and we used this list as a starting point.

From this list, we filtered all verbs (Step 1) and then filtered intransitive
verbs using LDOCE Online mentioned in Section 2.2 (Step 2). As shown in Table
2, limiting the analysis to intransitive verbs decreases the number of target
verbs, but it is a necessary procedure. This is because most English verbs
exhibit properties of both intransitive and transitive verbs [21].

Then, from the filtered intransitive verbs, we selected only those that
appeared more than 200 times in each of EEBO, COHA, and Google Books (Step 3).
This is the strictest condition for our quantification. We also examined a
mild condition (30 times in each dataset), as described later. We set these
values based on our pre-investigation results; 1) no target verbs remain when
the threshold is larger than 200, and therefore we set it as the upper bound;
2) when the threshold is below 30, the computed frequencies of less popular
verbs around the threshold are just unreliable and therefore it was set as the
lower bound.

Finally, we selected only those intransitive verbs with the be+PP frequency of
0.5 or greater in our oldest corpus EEBO. This eliminates intransitive verbs
for which the have+PP form was already dominant in earlier English. This
process aims to ensure that, for all the target verbs, the evolution of the
perfect began in the period covered by EEBO, COHA, and Google Books.

Table 2: Steps for the frequency-based verb selection. Step | Procedure | Selected verbs
---|---|---
1 | Select verbs from a full list of contemporary words in COCA | 5,764
2 |  | Select intransitive verbs based on LDOCE Online
---
719
3 |  | Select verbs from 2 that occurred more than 200 times
---
in the three datasets
46
4 |  Select verbs from 3 with more than 50% of be+PP in EEBO | 13

After these steps, we obtained 13 intransitive verbs that were in the be+PP
form in earlier English (EEBO) and sufficiently frequent in the three data
sources (1473–2000).

#### (B) Verb group selected based on prior studies

In addition to Group A, we investigated additional example verbs. Among
intransitive verbs that were excluded using the frequency-based method in verb
group A, verbs that were already considered in prior studies would be examined
separately from those of Group A. Among verbs listed in LDOCE Online with both
intransitive and transitive usage, prior studies [12, 14] indicated that the
six verbs were frequently used in the perfect with BE (see Group B in Table
3). Therefore, we analysed them separately from Group A and expected to find
evolutionary trajectories different from those of Group A, which may give us
additional insights into the cultural evolution of grammar.

Table 3: Target verbs. Group | Selected verbs
---|---
A |  | arrive, bound, come, creep, degenerate, expire, fall, insist, look,
---
rise, stay, tumble, vanish
B | ascend, become, depart, descend, escape, go

###  2.4 Detection of evolutionary forces

Unidirectional selection, in which new forms replace older forms over
generations, often (but not always) produces an S-shaped growth curve when
viewed as a change in allele frequency in a population [22], which is often
accepted as evidence of directional force favouring one variant over others.
The historical changes in the selection of auxiliary verbs be/have+PP can also
be expressed mathematically within the same framework. For detecting
evolutionary forces, there are two established methods, the Frequency
Increment Test (FIT) [8, 23] and a neural network-based classification (TSC)
[24]. Both methods are based on the Wright-Fisher model [25], which describes
the drift dynamics of two competing types in a population of fixed size $N$
and uses it as the null model.

The FIT is a method to detect selection in time series data, applying it to
population genetics experiments and historic DNA samples [23]. The FIT rejects
‘random drift’ (the null hypothesis) when the distribution of increments in
the frequency of normalised alleles shows a mean value that deviates
significantly from zero, which suggests the possibility of evolutionary forces
or ‘selection’. However, the FIT has several drawbacks as reported in the
literature [26]: it can be sensitive to how the corpus/dataset is organised
into temporal segments (i.e., binning), and it assumes the normality of data,
which the real data often violate.

To resolve these problems, Karsdorp et al. (2020) proposed a deep neural
network model, in which the problem of detecting evolutionary forces is
formulated as a binary classification task for a given time series data. This
is called the neural time series classification (TSC). This model was trained
on time series of cultural change simulated by the Wright-Fisher model. They
demonstrated that the neural TSC could resolve problems mentioned above
associated with the FIT: it has robustness for specific binning methods, and
the normality assumption of frequency increments does not play a role in this
model. Nevertheless, the neural TSC can consistently and accurately
distinguish time series produced by random drift from time series subject to
selection pressure.

Therefore, in this study, we used the neural TSC library [27] to detect
evolutionary forces underlying the transition from be+PP to have+PP. In
addition, for reference, we also used the FIT in our datasets. However, as
explained later (and in Supporting Material), our post-hoc power analysis
shows that the data size is insufficient to apply the FIT.

###  2.5 Data processing

Figure 1: Data processing. (A) Construction of the frequency time series for a
target verb using three data sources, in which the years represent the time
range covered by each corpus in our analysis (B) Binning method. The data size
per bin should be comparable.

Figure 1 shows schematic illustrations of our data processing for detecting
evolutionary forces. As explained in Section 2.3, we restricted our search
queries to the basic form of be/have+PP in EEBO, COHA, and Google Books to
prevent false positives (i.e., be+PP for the passive) and accurately compute
their frequencies of matches across time for each target verb. Since the
spellings of (auxiliary) verbs vary over history, we referred to the formats
found in EEBO and formulated lists (Tables S1–3 in Supplementary Material) to
perform a comprehensive search of be/have and the past participle of each
verb.

When merging the search results of be/havee+PP from the three datasets, it
should be kept in mind that, for EEBO and COHA, we can retrieve the patterns
that match the be/have+PP construction, along with the year in which the
sentence appeared (see Fig. 1A upper middle table). By aggregating these data,
we can construct a frequency time series given a target verb (raw counts of
occurrences per year). As mentioned before, with the Google Books Ngram
Viewer, we can only obtain the percentages of matches adjusted by year
(relative frequencies of occurrences per year) without the matched patterns.
Furthermore, as Google Books is more than three orders of magnitude larger
than those of the other two corpora (see Table 1), we must properly adjust the
search results by ‘scaling’.

Google Books and COHA have an overlapping period between 1810 and 2000. We,
therefore, took advantage of this fact to scale the search results from Google
Books so that the resulting frequencies roughly matched those from COHA.
Specifically, for the years 1810–2000, we focused on all verbs to find the
average frequency of be/have+PP in COHA ($f_{C}$), and then obtained the
average frequency in Google Books ($f_{G}$) to estimate a scaling constant of
$C:=f_{C}/f_{G}$. The multiplication of Google Books search results by $C$ was
used as data to complement the periods not covered by EEBO and COHA (i.e.,
1700–1810). More specifically, we reconstructed the search results for Google
Books in as those from EEBO and COHA (like the one in Fig. 1A upper middle)
based on the scaled frequency. No scaling was applied to EEBO as the corpus
size was comparable to COHA, and there was no overlapping period between these
corpora.

Figure 1B illustrates the binning for constructing relative frequency time
series of be/have+PP, which used the same setting reported in [8], where the
bin size was set to $\log N$ ($N$ being the total counts). Then, we split the
data so that each bin had approximately the same data size, which is a
necessary treatment because the neural TSC (and FIT) assumes the frequency
changes within the same population (i.e., the Write-Fisher model [25] as the
null model). The median of the year data in each bin was used as the time for
the bin.

The frequency time series of be/have+PP are shown in Fig. S1 and S2 in
Supplementary Material. For both Groups A and B, we can confirm that for all
verbs, the trend of frequency increase during the period of overlap between
Google Books and COHA (1810–2000) is consistent, and the frequency at the
border between EEBO and Google Books is mostly consistent. In other words,
this indicates that the Google Books results can properly complement the other
two corpora through scaling and binning. Although we are aware of little gaps
between 1700 and 1750 in the combined time series, this may not considerably
affect the results of the neural TSC, because the frequency time series were
smoothed by binning and the inflection point seems after these gaps, which was
about 1800 according to the literature [15, 16].

##  3 Results

Figure 2 shows the historical changes in the frequency of be/have+PP
constructions for the 13 verbs in Group A that were selected based on
frequency. For all verbs except bound, there was a clear increase in frequency
from be+PP to have+PP. While be+PP was dominant with most verbs before 1600,
there was a sharp increase in the frequency of have+PP between 1750 and 1800.
This result clearly shows that the form of have+PP became dominant during the
cultural evolution of grammar.

Furthermore, we performed the same analysis on the six intransitive verbs from
Group B. These verbs were analysed in prior work and excluded from the
selection criteria of verb Group A. As shown in Fig. 3, similar evolutionary
trends can be seen as in Group A; that is, be+PP transitions to have+PP. We
also note that the earlier onset (‘escape’) and the later onset (‘descend’)
are identified in the frequency increase, suggesting individual differences
between verbs.

Figure 2: Historical frequency changes of be/have+PP for the 13 verbs in Group A. Figure 3: Historical frequency change of be/have+PP for the six verbs in Group B. Table 4: Neural TSC result for Group A and B. Verb | Group | Probability
---|---|---
arrive | A | 1.00
come | A | 1.00
creep | A | 1.00
degenerate | A | 1.00
expire | A | 1.00
fall | A | 1.00
insist | A | 1.00
rise | A | 1.00
stay | A | 1.00
vanish | A | 1.00
ascend | B | 1.00
become | B | 1.00
depart | B | 1.00
escape | B | 1.00
look | A | 1.00
tumble | A | 1.00
descend | B | 0.78
bound | A | 0.07
go | B | 0.01

Were the frequency changes of be/have+PP due to random drift or directional
forces at play? The previous studies [15, 16] did not address this question
quantitatively. To examine this, we applied the neural TSC to all the verbs in
Group A and B. Table 4 summarizes the results, where if the probability is
above 0.5, it is deemed as ‘selection’. Of the 19 verbs, 17 verbs were
classified as ‘selection’. The verbs classified as random drift were ‘bound’
from Group A and ‘go’ from Group B. This is visually confirmed from Fig. 2 and
3; these verbs did not exhibit an apparent increasing trend or a tipping point
in the historical frequency changes of be/have+PP. It is likely that these
verbs were not classified as ‘selection’ because the increase in have+PP was
suppressed due to the presence of the be+PP passive usage (e.g., ‘England is
bounded on the north by Scotland.’ Is bounded here is passive.) or the
adjective usage of PP (e.g., ‘His money was gone and his health broken.’ Gone
here can be regarded as an adjective.). Note that ‘bound’ is categorised as an
intransitive verbs by LDOCE.

For reference, we also conducted the FIT for our datasets. The result is shown
in Table S4 in Supplementary material. Our post-hoc power analysis revealed
that the estimated power $d<0.8$ for all verbs, indicating that the data size
was insufficient for this statistical test [28], even though we used the
largest data sources. Therefore, we did not use the FIT results for
discussion.

##  4 Discussion

### Evolutionary forces in the cultural evolution of grammar

In modern English, the perfect takes the form of have+PP, and this grammatical
rule applies to most verbs. However, in earlier English, there were verbs
which formed the perfect with be+PP. As mentioned earlier, the prior studies
have addressed this phenomenon at the aggregated level; on average, the
frequency of the HAVE perfect started increasing in Late Middle English
(1420–1500), and that of the BE perfect started decreasing in the 19th century
(1810–1861) [15, 16]. Similar trends can be observed in Fig. 2 and 3.

In this study, we quantified evolutionary forces by applying a robust neural
network-based classification (the neural TSC) to the large-scale data sources:
EEBO, COHA, and Google Books with scaling. We found that most verbs in Group A
and B clearly exhibited an increase in the frequency of have+PP, among which
17 verbs were classified as ‘selection’ by the neural TSC. Two exceptional
cases are explained as the existence of the be+PP passive usage or the
adjective usage of PP, which can suppress the increase in the frequency of
have+PP. They could cause the neural TSC to fail in classification. Given the
fact that most verbs in Group A and B exhibited substantial frequency
increases in the form of have+PP, it is unlikely that the English perfect
evolved through random drift; selection might have played a major role in the
evolution of grammar. Although previous studies (described in Sec 1.2)
identified linguistic and extralinguistic factors accounting for the rise of
HAVE and the corresponding decline of BE with intransitives, our analysis
allows us to not only reconfirm the tipping point of the BE to HAVE transition
at an aggregated level but also provide evolutionary speeds, patterns, and
dynamics for each verb.

The role of random drift in language evolution has been emphasised in a
previous study [8], but our findings are contrasting in that respect.

We also examined the same experiment using a different setting. When lowering
the threshold for an expanded selection of target verbs to those appearing at
least 30 times in each dataset, comparable results were obtained. Most of the
target verbs (33/36) were classified as ‘selection’ and the exceptional verbs
(‘meddle’ in addition to ‘bound’ and ‘go’) are all explainable. For example,
consider the sentence, ‘This election was meddled with by the Russians.’ Here,
meddled is passive, although ‘meddle’ is classified as an intransitive verb by
LDOCE Online. Thus, this result further supports the conclusion that natural
selection was the major force for the transition from be+PP to have+PP.

### Limitations

There are limitations and several potential issues to be considered in our
study. First, our conclusion is based on the neural TSC that was trained on
time series of cultural change simulated by the Wright-Fisher model (as the
null model). Thus, if the frequency time series of a verb violated this
assumption, the neural TSC could not detect evolutionary forces even if it
existed, although the same applies to the other methods, including the FIT.
Similar to Newberry et al. (2017), the combined datasets of EEBO, COHA, and
Google Books include British and American English in a mixture, which might
affect the raw counts of matched search results because the pace at which
be/have+PP evolved could be different in British and American English, but the
overall tendency holds as we computed the relative frequency of be/have+PP. In
addition, there is a ‘bias’ in Google Books resulting from the increasing
inclusion of scientific texts [29]. Verbs that are commonly used in scientific
contexts differ from those that are prevalent in other contexts, and as a
result, the relative frequency of be+PP / have+PP may also vary by genre.
Therefore, it is crucial for future research to acknowledge the significance
of genre as a factor in the analysis.

Second, we showed that selection may have played a role in the evolution of
the English perfect in many verbs. However, the linguistic reasons for the
rise of have+PP and the fall of be+PP are excluded from our consideration,
although linguistic factors, as well as extralinguistic ones, are relevant to
the dissipation of be+PP, as explained in the previous studies. Moreover, it
is known that throughout its evolution, the perfect has become capable of
expressing not only ‘completion’ and ‘result’, but also ‘experience’ and
‘continuity’. We speculate that the development of such new functions might be
related to the shift of auxiliary verbs from BE to HAVE. This may be because
functional differentiation between BE and HAVE could be adaptive to avoid
expressive confusion. That is, BE is often used for the passive rather than
the perfect, and therefore be+PP has always been at risk of being mistaken for
the passive. The evolution of the have+PP rule reduces such a risk, allowing
for the expression of complex meanings. This is still at the stage of a
hypothesis and needs further investigation. Going forward, we should conduct
further studies on the relation between the evolution of the perfect and the
functional developments within the categories of Aspect and Modality.
Moreover, with the same analytical framework, languages other than English
should be analysed to investigate the universal dynamics of linguistic
diversity resulting from cultural evolution.

### Summary

In summary, our study confirmed the findings of prior work [12, 13, 14, 15,
16] using the larger data sources, and further provided evidence of the
possibility of evolutionary forces at work in the transition from be+PP to
have+PP for each verb. Whether the evolution of the perfect is due to
selection or drift in other languages remains an important future research
question. Our study provides detailed descriptions of processes helpful for
conducting such attempts in multiple languages, which is intrinsic to
developing the theoretical framework of the cultural evolution of grammar.

## Data Accessibility

The data sources used in this study are available:
EEBO (https://www.english-corpora.org/eebo/);
COHA (https://www.english-corpora.org/coha/);
Google Books (https://books.google.com/ngrams/).
The code and process data are available at https://bit.ly/3AVP3F9 and
http://bit.ly/3TMQrDy. (This will be moved to a public repository on
publication.)

## Author Contributions

S.O.: Data curation, formal analysis, investigation, methodology, validation,
visualization, writing—review and editing. M.H.: conceptualization, funding
acquisition, investigation, writing—review and editing; K.S.:
conceptualization, funding acquisition, investigation, writing—original draft,
writing—review and editing. MH and KS designed the study. SO and KS analysed
data. All the authors have written and revised the manuscript.

## Competing interests

We declare that we have no competing interests.

## Ethical statements

This article does not contain any studies involving human participants.

## Funding

This work was supported by MEXT/JSPS KAKENHI Grant Numbers #4903, JP17H06379
and JP17H06383.

## Acknowledgements

We would like to thank the members of the Evolinguistic project for their
fruitful discussions.

## References

  * [1] Dediu, D. _et al._ Cultural evolution of language.  In _Cultural evolution: Society, technology, language, and religion. Strüngmann Forum Reports, vol. 12_ , 303–332 (MIT Press, 2013).
  * [2] Mesoudi, A.  _Cultural evolution_ (University of Chicago Press, 2011).
  * [3] Hauser, M. D., Chomsky, N. & Fitch, W. T.  The Faculty of Language: What Is It, Who Has It, and How Did It Evolve?  _Science_ 298, 1569–1579 (2002).
  * [4] Kirby, S., Cornish, H. & Smith, K.  Cumulative cultural evolution in the laboratory: An experimental approach to the origins of structure in human language.  _Proceedings of the National Academy of Sciences_ 105, 10681–10686 (2008).
  * [5] Blythe, R. A.  Neutral evolution: A null model for language dynamics.  _Advances in Complex Systems_ 15, 1–20 (2012).
  * [6] Kirby, S., Griffiths, T. & Smith, K.  Iterated learning and the evolution of language.  _Current Opinion in Neurobiology_ 28, 108–114 (2014).
  * [7] Lieberman, E., Michel, J.-B., Jackson, J., Tang, T. & Nowak, M. A.  Quantifying the evolutionary dynamics of language.  _Nature_ 449, 713–716 (2007).
  * [8] Newberry, M. G., Ahern, C. A., Clark, R. & Plotkin, J. B.  Detecting evolutionary forces in language change.  _Nature_ 551, 223–226 (2017).
  * [9] Aranovich, R.  _Split Auxiliary Systems: A Cross-linguistic Perspective_.  Typological studies in language (John Benjamins Pub., 2007).
  * [10] Sorace, A.  Gradients in auxiliary selection with intransitive verbs.  _Language_ 76, 859–890 (2000).
  * [11] Ackema, P. & Sorace, A.  Auxiliary selection.  In _The Wiley Blackwell Companion to Syntax, Second Edition_ , 1–32 (John Wiley & Sons, Ltd, 2017).
  * [12] Kytö, M.  Be/have + past participle: The choice of the auxiliary with intransitives from Late Middle to Modern English.  In Rissanen, M., Kytö, M. & Heikkonen, K. (eds.) _Corpus-based Studies in Linguistic Variation and Genre Styles_ , 17–86 (De Gruyter Mouton, Berlin, Boston, 2011).
  * [13] Huber, J.  Counterfactuality and aktionsart: Predictors for BE vs. HAVE + past participle in Middle English.  In Claridge, C. & Bös, B. (eds.) _Developments in English historical morpho-syntax_ , 149–173 (John Benjamins, Amsterdam, 2019).
  * [14] Rydén, M. & Brorström, S.  _The be/have variation with intransitives in English. With special reference to the Late Modern Period_ , vol. 70 of _Stockholm studies in English_ (Almqvist och Wiksell, 1987).
  * [15] McFadden, T. & Alexiadou, A.  Perfects, Resultatives, and Auxiliaries in Earlier English.  _Linguistic Inquiry_ 41, 389–425 (2010).
  * [16] McFadden, T.  On the disappearance of the BE perfect in Late Modern English.  _Acta Linguistica Hafniensia_ 49, 159–175 (2017).
  * [17] Lodge, T.  Early English Books Online.  https://www.english-corpora.org/eebo/ (2017).
  * [18] Davies, M.  Corpus of Historical American English.  https://www.english-corpora.org/coha/ (2010).
  * [19] Google.  Google Ngram Viewer.  http://books.google.com/ngrams/datasets (2012).
  * [20] Education, P.  Longman Dictionary of Contemporary English Online.  https://www.ldoceonline.com/jp/dictionary/, https://www.ldoceonline.com/jp/dictionary/english-japanese/.
  * [21] van Gelderen, E.  _The Diachrony of Verb Meaning:Aspect and Argument Structure_ (Routledge, 2018).
  * [22] Blythe, R. A. & Croft, W.  S-curves and the mechanisms of propagation in language change.  _Language_ 88, 269–304 (2012).
  * [23] Feder, A. F., Kryazhimskiy, S. & Plotkin, J. B.  Identifying Signatures of Selection in Genetic Time Series.  _Genetics_ 196, 509–522 (2014).
  * [24] Karsdorp, F., Manjavacas, E., Fonteyn, L. & Kestemont, M.  Classifying evolutionary forces in language change using neural networks.  _Evolutionary Human Sciences_ 2, 1–16 (2020).
  * [25] Ewens, W. J.  _Mathematical Population Genetics 1: Theoretical Introduction_.  Interdisciplinary Applied Mathematics (Springer, New York, 2012).
  * [26] Karjus, A., Blythe, R. A., Kirby, S. & Smith, K.  Challenges in detecting evolutionary forces in language change using diachronic corpora.  _Glossa: a journal of general linguistics_ 5, 1–25 (2020).
  * [27] Folgert Karsdorp.  nnfit.  https://github.com/fbkarsdorp/nnfit.
  * [28] Cohen, J.  _Statistical Power Analysis for the Behavioral Sciences_ (Taylor & Francis, 2013).
  * [29] Pechenick, E. A., Danforth, C. M. & Dodds, P. S.  Characterizing the Google Books Corpus: Strong Limits to Inferences of Socio-Cultural and Linguistic Evolution.  _PLOS ONE_ 10, 1–24 (2015).