# Spectral function for 4He using the Chebyshev expansion in coupled-cluster
theory

J. E. Sobczyk  Institut für Kernphysik and PRISMA+ Cluster of Excellence,
Johannes Gutenberg-Universität, 55128 Mainz, Germany  S. Bacca  Institut für
Kernphysik and PRISMA+ Cluster of Excellence, Johannes Gutenberg-Universität,
55128 Mainz, Germany  Helmholtz-Institut Mainz, Johannes Gutenberg-Universität
Mainz, D-55099 Mainz, Germany  G. Hagen  Physics Division, Oak Ridge National
Laboratory, Oak Ridge, TN 37831, USA  Department of Physics and Astronomy,
University of Tennessee, Knoxville, TN 37996, USA  T. Papenbrock  Department
of Physics and Astronomy, University of Tennessee, Knoxville, TN 37996, USA
Physics Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA

###### Abstract

We compute spectral function for 4He by combining coupled-cluster theory with
an expansion of integral transforms into Chebyshev polynomials. Our method
allows to estimate the uncertainty of spectral reconstruction. The properties
of the Chebyshev polynomials make the procedure numerically stable and
considerably lower in memory usage than the typically employed Lanczos
algorithm. We benchmark our predictions with other calculations in the
literature and with electron scattering data in the quasi-elastic peak. The
spectral function formalism allows one to extend ab-initio lepton-nucleus
cross sections into the relativistic regime. This makes it a promising tool
for modeling this process at higher energy transfers. The results we present
open the door for studies of heavier nuclei, important for the neutrino
oscillation programs.

##  I Introduction

Lepton-nucleus cross sections are not only a invaluable tool to investigate
the nuclear dynamics with clean electroweak probes Bacca and Pastore (2014),
but have also become a hot topic in the short- and long-baseline neutrino
programs aiming at extracting neutrino oscillation parameters Alvarez-Ruso
_et al._ (2018); Ankowski  _et al._ (2022); Ruso  _et al._ (2022). Recently,
we initiated a theory program that addresses low- and intermediate- energy
lepton-nucleus scattering from first principles by combining the Lorentz-
integral-transform with coupled-cluster theory (LIT-CC) to compute many-body
response functions Bacca  _et al._ (2013) and lepton-nucleus cross sections.
In this approach, the final state interaction is described consistently with
the initial state interaction by the same Hamiltonian rooted into quantum
chromodynamics, see Refs. Acharya and Bacca (2020); Sobczyk  _et al._ (2020,
2021). An example is the description of the longitudinal quasi-elastic peak of
40Ca, see Ref. Sobczyk  _et al._ (2021), and this paves the way for further
investigations. Despite its success, this approach also has its limitations:
the formalism is based on the non-relativistic theory and at the moment is
capable of predicting only inclusive cross-sections. Efforts to include the
spectrum of the outgoing nucleons below and above pion production are still
lacking in ab-initio calculations. Several other approximations and
phenomenological methods instead offer a way to answer such questions Nieves
_et al._ (2011); Martinez  _et al._ (2006); Amaro  _et al._ (2005); Leitner
_et al._ (2009); Benhar  _et al._ (2005), chief among them being the spectral
functions formalism. While spectral functions can be computed
phenomenologically Benhar  _et al._ (1994); Nieves and Sobczyk (2017); Buss
_et al._ (2007), calculations were performed recently within the self-
consistent Green function (SCGF) method Rocco and Barbieri (2018); Barbieri
_et al._ (2019) using a similar Hamiltonian as in Ref. Sobczyk  _et al._
(2021) for the initial state. In the past also the LIT method combined with
hyperspherical harmonics was used to obtain the proton spectral function of
4He Efros  _et al._ (1998).

The main advantage of the spectral function formalism lies in the possibility
of detaching the high-energy physics from the ground-state properties of the
nucleus, under the assumption that the final state interactions can be
neglected. This not only allows one to make predictions for the quasi-elastic
peak using the relativistic kinematics and currents, but this approach can
also be used at higher energies, e.g., above the pion production threshold.
Hence, developing efficient ab-initio methods to compute spectral functions
deserves attention, which goes beyond the mere fact that their calculation in
a many-body system is per se an interesting and challenging task. In this
work, we present an approach to the computation of spectral functions which
opens up the possibility of using ab-initio many-body methods in the high-
energy regime.

The reconstruction of nuclear response functions requires information about
the excited states of the system, but usually these are not easily accessible.
To circumvent this issue, the problem has been often reformulated by computing
integral transforms of the response function, with Lorentz and Laplace kernels
being popular choices Efros  _et al._ (1994, 2007); Carlson and Schiavilla
(1992); Carlson  _et al._ (2015); Lovato  _et al._ (2016). The computation of
the integral transform requires one to (only) solve a bound-state problem. The
inversion of the transform, needed to obtain the response function, has to be
performed numerically. While accurate results have been obtained for a variety
of electroweak observables Bacca  _et al._ (2013, 2014); Lovato  _et al._
(2016, 2020); Sobczyk  _et al._ (2021), the inversion introduces an additional
numerical error and is most stable when the response function exhibits only
one or two broad peaks Bacca  _et al._ (2002, 2004). This scenario was
recently explored using machine-learning techniques Raghavan  _et al._ (2021).
Spectral functions often have a more complicated structure and this makes the
inversion of the integral transform difficult. Therefore, we propose here to
use a different approach that is based on the Chebyshev expansion of the
integral kernel (ChEK) introduced in Refs. Roggero (2020); Sobczyk and Roggero
(2021). Although it relies on the idea of the integral transform, it does not
require its inversion. Moreover, for a given desired resolution of the
spectral reconstruction it allows one to estimate an uncertainty.

The information about the (discrete) spectrum of excited states in a many-body
system can be retrieved in various ways. The nuclear theory community is
familiar with the Lanczos orthogonalization procedure Lanczos (1950), which
for example is used in the LIT-CC method. An alternative approach, developed
in the field of condensed-matter physics, is the kernel polynomial method
(KPM) Weiße  _et al._ (2006). As for the KPM, the approach of this paper is
also based on an expansion in Chebyshev polynomials.

This paper is organized as follows. In Sec. II, we review how the lepton-
nucleus interaction in the quasi-elastic peak can be expressed in terms of
spectral functions using the impulse approximation. In Sec. III, we present
the theoretical framework for our calculations of spectral functions. We
validate our method in 4He, paying special attention to the center-of-mass
problem in Sec. IV, and finally we conclude in Sec. V.

##  II Electron-nucleus scattering

Let us consider the process

$e(k)+A(p_{0})\to e^{\prime}(k^{\prime})+f(p_{f})\,,$ (1)

where an incoming electron with four-momentum $k=(E_{k},\mathbf{k})$ is
scattered off a nucleus $A$, producing an outgoing electron with four-momentum
$k^{\prime}=(E_{k}^{\prime},\mathbf{k}^{\prime})$ and a final (in general
multi-particle) state $f$. The four-momentum transfer is
$q\equiv(\omega,\mathbf{q})=k-k^{\prime}$. In the Born approximation, the
electron interacts via the exchange of a single $\gamma$.

The inclusive cross section of this process can be written in terms of
leptonic and hadronic tensors as

$\frac{d\sigma}{d\omega
d\cos\theta}=\left(\frac{\alpha}{q^{2}}\right)^{2}\frac{|\mathbf{k}|}{|\mathbf{k}^{\prime}|}L_{\mu\nu}W^{\mu\nu}\,,$
(2)

with the angle of the outgoing lepton being $\theta$, $\alpha\approx 1/137$,
is the fine structure constant. The lepton tensor is

$L_{\mu\nu}=2[k_{\mu}k_{\nu}^{\prime}+k_{\mu}^{\prime}k_{\nu}-g_{\mu\nu}(kk^{\prime})]\,.$
(3)

The nuclear structure information is encoded in the hadron tensor

$W^{\mu\nu}=\sum_{f}\delta^{4}(p_{0}+q-p_{f})\langle
0|\left(J^{\mu}\right)^{\dagger}|\Psi_{f}\rangle\langle\Psi_{f}|J^{\nu}|0\rangle\,,$
(4)

where the current $J^{\mu}$ corresponds to the electromagnetic process.

In the following, we will focus on the quasi-elastic mechanism, for which the
interaction takes place on a single nucleon, kicking off a nucleon from the
remaining $(A-1)$ nucleus in the final state. The electromagnetic current is a
sum of one-body contributions which in the second quantization form is given
by

$J^{\mu}=\sum_{\alpha,\beta}\langle\beta|j^{\mu}|\alpha\rangle
a_{\beta}^{\dagger}a_{\alpha}\,.$ (5)

Here, the initial and final nucleon states are labelled by $\alpha$ and
$\beta$, respectively. Within the spectral function formalism, we will use the
fully relativistic current in the matrix element treating the initial and
final nucleons as free states

$\langle p+q|j^{\mu}|p\rangle=\bar{u}(p+q)V^{\mu}u(p)\ .$ (6)

The current $j^{\mu}$ has a vector structure and $u$ denotes a Dirac spinor.
Constructing the most general form of $V^{\mu}$ using the available four-
vectors, we have

$V^{\mu}=F_{1}\gamma^{\mu}+\frac{F_{2}}{2m}i\sigma^{\mu\nu}q_{\nu}\,.$ (7)

We use $F_{1}^{n,p}$, $F_{2}^{n,p}$ parametrized as in Ref. Bradford  _et al._
(2006).

###  II.1 The Impulse Approximation

At relatively large momentum transfer $\mathbf{q}$, one can assume that the
struck nucleon is decoupled from the nuclear $(A-1)$ system, i.e., that the
final-state interaction can be neglected. Within this impulse approximation,
the final nuclear state factorizes as

$|\Psi_{f}\rangle\longrightarrow
a_{p^{\prime}}^{\dagger}|\Psi_{A-1}\rangle\,,$ (8)

where a plane-wave state $a_{p^{\prime}}$ with momentum $\mathbf{p^{\prime}}$
and energy $E_{p^{\prime}}$ is added on top of the final $(A-1)$ system 111For
simplicity we suppress spin and isospin indices.. Using the current of Eq.
(5), the one-body matrix element can be factorized as

$\displaystyle\langle\Psi_{f}|J^{\mu}|0\rangle$
$\displaystyle\to\sum_{\alpha,\beta}\langle\beta|j^{\mu}|\alpha\rangle\langle\Psi_{A-1}|a_{p^{\prime}}a_{\beta}^{\dagger}a_{\alpha}|0\rangle$
$\displaystyle=\sum_{\alpha,\beta}\langle\beta|j^{\mu}|\alpha\rangle\langle\Psi_{A-1}|\delta_{\beta\,p^{\prime}}a_{\alpha}-a_{\beta}^{\dagger}a_{p^{\prime}}a_{\alpha}|0\rangle$
$\displaystyle\approx\sum_{\alpha}\langle
p^{\prime}|j^{\mu}|\alpha\rangle\langle\Psi_{A-1}|a_{\alpha}|0\rangle$
$\displaystyle=\int\frac{d^{3}\mathbf{p}}{(2\pi)^{3}}\langle
p^{\prime}|j^{\mu}|p\rangle\sum_{\alpha}\langle\mathbf{p}|\alpha\rangle\langle\Psi_{A-1}|a_{\alpha}|0\rangle\,,$
(9)

where the approximation in the third line assumes that the struck nucleon at
the interaction vertex is exactly the one which is ejected from the nucleus
Dickhoff and Van Neck (2005) and in the last line we insert a complete set of
states $\int d^{3}p/(2\pi)^{3}|p\rangle\langle p|$.
$\langle\mathbf{p}|\alpha\rangle$ are single-particle wavefunctions in
momentum space. The process is shown schematically in Fig. 1.

Figure 1: Impulse approximation: the electroweak interaction takes place on a
single nucleon which afterwards does not interact with the spectator
$|\Psi_{A-1}\rangle$ system.

The recoil energy $E_{f}^{kin}$ of the $(A-1)$ system is negligible for heavy
nuclei, and the excitation is given by the energy conservation

$E_{f}^{A-1}=\omega+E_{0}-E_{p^{\prime}}-E_{f}^{kin}$ (10)

with the initial-state energy $E_{0}$.

Substituting Eq. (II.1) in Eq. (4), the incoherent contribution to the hadron
tensor becomes

$\displaystyle W^{\mu\nu}(q)=$
$\displaystyle\int\frac{d^{3}\mathbf{p}d^{3}\mathbf{\tilde{p}}}{(2\pi)^{6}}\sum_{\alpha,\alpha^{\prime}}\sumint_{\Psi_{A-1}}\,\langle
p|\left(j^{\mu}\right)^{\dagger}|p^{\prime}\rangle\langle
p^{\prime}|j^{\nu}|\tilde{p}\rangle$
$\displaystyle\times\langle\mathbf{p}|\alpha\rangle^{\dagger}\langle\mathbf{\tilde{p}}|\alpha^{\prime}\rangle\langle
0|a_{\alpha}^{\dagger}|\Psi_{A-1}\rangle\langle\Psi_{A-1}|a_{\alpha^{\prime}}|0\rangle$
$\displaystyle\delta(\omega-
E_{p^{\prime}}-E_{f}^{A-1}-E_{f}^{kin}+E_{0})\,\,.$ (11)

From the momentum conservation at the single nucleon vertex
$\mathbf{p}=\mathbf{\tilde{p}}=\mathbf{p^{\prime}}-\mathbf{q}$. Furthermore,
the spin state of $a_{\alpha}$ and $a_{\alpha^{\prime}}$ coincide due to
charge conservation and the assumption that the nuclear ground state has spin
zero. Finally, the last step of the factorization separates the energy
conservation at the vertex from the excitation of the residual system by
introducing the energy $E$ needed to remove a nucleon with momentum
$\mathbf{p}$ from the ground state as

$\displaystyle\delta(\omega+E_{0}-E_{p^{\prime}}-E_{f}^{A-1}-E_{f}^{kin})$
$\displaystyle=\int
dE\delta(\omega+E-E_{p^{\prime}}-E_{f}^{kin})\delta(E+E_{f}^{A-1}-E_{0})\,.$
(12)

Using this equation and introducing explicitly the isospin dependence, the
hadron tensor is

$\displaystyle W^{\mu\nu}(q)$
$\displaystyle=\int\frac{d^{3}\mathbf{p}}{(2\pi)^{3}}dE\big{[}S^{n}(\mathbf{p},E)w_{n}^{\mu\nu}(p,q)+$
$\displaystyle
S^{p}(\mathbf{p},E)w_{p}^{\mu\nu}(p,q)\big{]}\delta(\omega+E-E_{p+q}-E_{f}^{kin})\,,$
(13)

with $w_{n,p}^{\mu\nu}(p,q)=\langle p+q|j^{\mu}|p\rangle^{\dagger}\langle
p+q|j^{\nu}|p\rangle$ depending on the isospin of $|p\rangle$, where we have
introduced the hole spectral function

$\displaystyle S^{n,p}(\mathbf{p},E)$
$\displaystyle=\sum_{\alpha,\alpha^{\prime}}\sumint_{\Psi_{A-1}}|\langle
0|a_{\alpha}^{\dagger}|\Psi_{A-1}\rangle\langle\Psi_{A-1}|a_{\alpha^{\prime}}|0\rangle$
$\displaystyle\langle\mathbf{p}|\alpha\rangle^{\dagger}\langle\mathbf{p}|\alpha^{\prime}\rangle\delta(E+E_{f}^{A-1}-E_{0})\,.$
(14)

The spectral function gives the probability distribution of removing a nucleon
with momentum ${\mathbf{p}}$ from the target nucleus, leaving the residual
$(A-1)$ system with an energy $E_{0}-E$. For closed-shell nuclei, such as the
4He considered in this work, the spectral functions of spin-up and spin-down
nucleons coincide. We normalize spectral functions as

$\int\frac{d^{3}p}{(2\pi)^{3}}dES^{n(p)}(\mathbf{p},E)=N(Z)\,.$ (15)

In the relativistic regimes, the factors $m/E_{p}$ and $m/E_{p+q}$ should be
included to account for the implicit covariant normalization of the four-
spinors of the nucleons in the matrix elements of the current $j^{\mu}$. Hence
the hadron tensor finally becomes

$\displaystyle W^{\mu\nu}(q)$
$\displaystyle=\int\frac{d^{3}\mathbf{p}}{(2\pi)^{3}}dE\frac{m}{E_{p}}\frac{m}{E_{p+q}}$
$\displaystyle\big{[}S^{n}(\mathbf{p},E)w_{n}^{\mu\nu}(p,q)+S^{p}(\mathbf{p},E)w_{p}^{\mu\nu}(p,q)\big{]}$
$\displaystyle\times\delta(\omega+E-E_{p+q}-E_{f}^{kin})\,,$ (16)

where one can see that it can be calculated starting from the spectral
function for $(n)$ neutrons and $(p)$ protons.

We performed the factorization of the relativistic currents and the nuclear
ground state governed by non-relativistic dynamics. This way we can address
the processes occurring at high energy-momentum transfers. This procedure
introduces, however, some model dependence since we do not treat the
wavefunction and the currents on an equal footing. In a consistent description
we can either have the picture of a simple current interacting with a
complicated nucleus, or an alternative picture (and a continuum of approaches
between these two extremes) where the nucleus is simple (e.g., a product
state), and the current is complicated and consists of one- and two-body
terms. Very recently, authors of Ref. Tropiano  _et al._ (2021) presented a
detailed discussion of this subject. In particular they analyse how the high-
momentum behaviour of the wavefunctions depends on the resolution of employed
nuclear Hamiltonian. Their results can be applied to the momentum distribution
of the spectral functions. We leave the analysis of this effect, as well the
role of two-body currents within the factorization scheme, for the future
studies.

##  III Formalism

###  III.1 Green’s function and spectral function

The spectral functions Eq. (14) are defined through the imaginary part of a
propagator in a many-body system. Presently, we will consider only the hole
propagation of a state with quantum numbers $\alpha$ to a state with quantum
numbers $\beta$

$\begin{split}&G_{h}(\alpha,\beta,E)=\langle
0|a_{\beta}^{\dagger}\frac{1}{E-(E_{0}-\hat{H})-i\epsilon}a_{\alpha}|0\rangle\,,\\\
&\mathrm{Im}G_{h}(\alpha,\beta,E)=\\\ &-\pi\sumint_{\Psi_{A-1}}\langle
0|a_{\beta}^{\dagger}|\Psi_{A-1}\rangle\langle\Psi_{A-1}|a_{\alpha}|0\rangle\delta\big{(}E-(E_{0}-E_{\Psi})\big{)}\,.\end{split}$
(17)

The spectral functions can be retrieved from the imaginary part of Green’s
function summing over all the appropriate quantum numbers

$S(\mathbf{p},E)=-\frac{1}{\pi}\sum_{\alpha,\beta}\langle\mathbf{p}|\alpha\rangle\langle\mathbf{p}|\beta\rangle^{\dagger}\mathrm{Im}G_{h}(\alpha,\beta,E)\,,$
(18)

The reconstruction of $\mathrm{Im}G_{h}(\alpha,\beta,E)$ from Eq. (17)
requires a summation over all the excited states $|\Psi_{A-1}\rangle$, which
contains not only bound states but also continuum states. Within ab-initio
methods the calculated spectrum is typically discretized because of the
truncation of the many-body space. Continuum effects can be included via
complex scaling techniques Suzuki  _et al._ (2005); Carbonell  _et al._ (2014)
or in the Berggren basis. The latter idea has been recently applied to obtain
the microscopical optical potential from the coupled-cluster theory Rotureau
_et al._ (2017, 2018). These techniques are usually combined with the Lanczos
algorithm to construct tri-diagonal forms of large matrices, and thereby give
access to the extreme eigenvalues of the problem. Here, however, we will use
another approach described in the next section.

###  III.2 Chebyshev expansion of integral transform

Within the ChEK method we rephrase our problem: instead of reconstructing the
response we want to estimate observables which are expressed as the energy
integrals of the response. The method can be used in a general situation

$\Phi=\int d\omega f(\omega)R(\omega)\,,$ (19)

where $f(\omega)$ is any bound function defining the observable and
$R(\omega)$ is a response function – which in our case corresponds to
$\mathrm{Im}G_{h}(\alpha,\beta,\omega)$. Our strategy to approximate the
quantity in Eq. (19) consists in applying the integral transform
$\tilde{R}(E)$

$\tilde{\Phi}=\int dEf(E)\tilde{R}(E)\,,$ (20)

in such a way that we control the approximation error $|\Phi-\tilde{\Phi}|$.
Let us also notice, that the reconstruction of $\Phi$ does not require the
inversion of integral transform. In our case $\tilde{R}(E)$ is given by

$\begin{split}\mathrm{Im}\tilde{G}_{h}(\alpha,\beta,E)&\\\ &=\int
d\omega\mathrm{Im}G_{h}(\alpha,\beta,\omega)K(\omega,E)\\\
&=-\pi\sumint_{\Psi_{A-1}}\langle
0|a_{\beta}^{\dagger}|K\big{(}E_{\Psi},E-E_{0}\big{)}|a_{\alpha}|0\rangle\\\
&=-\pi\langle
0|a_{\beta}^{\dagger}|K\big{(}\hat{H},E-E_{0}\big{)}|a_{\alpha}|0\rangle\,.\end{split}$
(21)

The kernel $K(\omega,E)$ can be realized by various functions. Here, we will
apply the Gaussian kernel

$K(\omega,E;\lambda)=\frac{1}{\sqrt{2\pi}\lambda}\exp\left(-\frac{(\omega-E)^{2}}{2\lambda^{2}}\right)\,.$
(22)

Following Ref. Roggero (2020), we characterize the kernel as $\Sigma$-accurate
with $\Lambda$-resolution:

$\sup_{\omega_{0}\in[-1,1]}\sumint_{\omega_{0}-\Lambda}^{\omega_{0}+\Lambda}K(\omega_{0},E)dE\geq
1-\Sigma\,.$ (23)

With these definitions we can provide the uncertainty bound for
$|\Phi-\tilde{\Phi}|$, which will depend on the properties of the function $f$
and the kernel $K$.

Next we expand Eq. (21) into Chebyshev polynomials and truncate the number of
terms $N$. This truncation will introduce an additional error, as will be
explained later. The truncated kernel

$K(\omega,E)=\sum_{k=0}^{N}c_{k}(E)\mathcal{T}_{k}(\omega)\,,$ (24)

is expressed in terms of $\mathcal{T}_{k}(\omega)=\cos[k\arccos(\omega)]$
which follow a recursive relation

$\begin{split}&\mathcal{T}_{0}(x)=1;\,\,\,\,\,\mathcal{T}_{-1}(x)=\mathcal{T}_{1}(x)=x;\\\
&\mathcal{T}_{k+1}(x)=2x\mathcal{T}_{k}(x)-\mathcal{T}_{k-1}(x)\,.\end{split}$
(25)

Let us assume that the Hamiltonian norm is known, and that we are able to
rescale our problem $[E_{min},E_{max}]\rightarrow[-1,1]$. This allows us to
use Chebyshev polynomials (which are defined on the interval $[-1,1]$). The
Hamiltonian spectrum can be obtained, e.g., via the Lanczos algorithm. The
rescaling is given then by:

$\begin{split}&a=(E_{max}-E_{min})/2\,,\ \ \ \ \ b=(E_{max}+E_{min})/2\\\
&H:=(H-b)/a\,.\end{split}$ (26)

Combining Eqs. (21) and (24) we obtain

$\begin{split}\mathrm{Im}\tilde{G}_{h}(\alpha,\beta,E)&\\\
=&-\pi\sum_{k=0}^{N}c_{k}(E)\langle
0|a_{\beta}^{\dagger}\mathcal{T}_{k}\big{(}\hat{H}\big{)}a_{\alpha}|0\rangle\\\
\equiv&-\pi\sum_{k=0}^{N}c_{k}(E)\mu_{k}\,.\end{split}$ (27)

For simplicity we will abuse the notation and understand that
$\mathrm{Im}\tilde{G}_{h}(\alpha,\beta,E)$ has an implicit $N$ dependence.
Furthermore, the moments $\mu_{k}$ have an implicit dependence on $\alpha$ and
$\beta$. The moments of the expansion $\mu_{k}$ can be retrieved from a many-
body calculation, using the recursive relation from Eq. (25)

$\begin{split}&\langle\tilde{\Psi}_{0}|\equiv\langle
0|a_{\beta}^{\dagger}|\,,\ \ \ \ \ |\Psi_{0}\rangle\equiv
a_{\alpha}|0\rangle\,,\\\
&\langle\tilde{\Psi}_{k}|\equiv\langle\tilde{\Psi}_{k-1}|\hat{H}\ \ \ \ \
|\Psi_{k}\rangle=\hat{H}|\Psi_{k-1}\rangle\\\
&\mu_{0}=\langle\tilde{\Psi}_{0}|\Psi_{0}\rangle\,,\ \ \ \ \
\mu_{1}=\langle\tilde{\Psi}_{0}|\Psi_{1}\rangle\equiv\langle\tilde{\Psi}_{1}|\Psi_{0}\rangle\\\
&\mu_{k+1}=2\langle\tilde{\Psi}_{0}|\Psi_{k+1}\rangle-\mu_{k-1}\equiv
2\langle\tilde{\Psi}_{k+1}|\Psi_{0}\rangle-\mu_{k-1}\,.\end{split}$ (28)

In the $(k+1)$-th step only the $|\Psi_{k}\rangle$ (or
$\langle\tilde{\Psi}_{k}|$) state has to be known from the previous iteration.
Similarly to the Lanczos procedure, we iterate the action of Hamiltonian
$\hat{H}$. Here, however, no orthogonality restoration is needed at each step,
which makes the procedure faster and requires less memory. The coefficients
$c_{k}$ from Eq. (27) depend on the chosen kernel and their form can be found
in Ref. Sobczyk and Roggero (2021).

In the present case, we will define the function $f(\omega)$ in Eq. (19) as a
histogram bin centered at $\eta$ and a half-width $\Delta$

$f(\omega)\equiv
h_{\Delta}(\eta,\omega)=\bigg{\\{}\begin{matrix}0&|\eta-\omega|>\Delta\\\
1&\text{otherwise}\end{matrix}\;.$ (29)

We are then interested in approximating the histogram

$\begin{split}\mathrm{Im}G_{h}(\alpha,\beta;\eta,\Delta)=\int
h_{\Delta}(\eta,E)\mathrm{Im}G_{h}(\alpha,\beta,E)dE\end{split}$ (30)

using its integral transform (given in Eq. (27)) with a finite number of
Chebyshev moments $N$

$\begin{split}\mathrm{Im}\tilde{G}_{h}(\alpha,\beta;\eta,\Delta)=-\pi\sum_{k=0}^{N}\mu_{k}\int
h_{\Delta}(\eta,E)c_{k}(E)dE\,.\end{split}$ (31)

As shown in Ref. Sobczyk and Roggero (2021), we get

$\begin{split}\mathrm{Im}\tilde{G}_{h}(\Delta-&\Lambda)-\Sigma-2\gamma(\Delta-\Lambda)\\\
&\leq\mathrm{Im}G_{h}(\Delta)\\\
\leq\mathrm{Im}\tilde{G}_{h}&(\Delta+\Lambda)+\Sigma+2\gamma(\Delta+\Lambda)\,,\end{split}$
(32)

where we used a shorter notation
$\mathrm{Im}G_{h}(\alpha,\beta;\eta,\Delta)\equiv\mathrm{Im}G_{h}(\Delta)$.
The truncation error $\gamma$ depends on the number of moments $N$ and the
properties of the kernel

$\gamma=\sum_{k=N}^{\infty}c_{k}(E)\mathcal{T}_{k}(\omega)\,.$ (33)

The analytical expression for the bounds on $\gamma$ can be found in Eqs.
(B5), (B22) of Ref. Sobczyk and Roggero (2021) for the Lorentzian and Gaussian
kernels, respectively. As has been advocated in Ref. Sobczyk and Roggero
(2021), the Gaussian integral transform has better convergence properties and
will therefore be used in the present calculations.

Eq. (32) is the master equation which will ultimately allow us to reconstruct
the spectral functions defined in Eq. (18) as a histogram. It gives the error
bound for $\mathrm{Im}G_{h}(\Delta)$ depending on the characteristics of the
kernel, $\Lambda$ and $\Sigma$, and the number of Chebyshev moments $N$ (which
enter both $\gamma$ and $\mathrm{Im}\tilde{G}_{h}(\Delta\pm\Lambda)$).

It is important to notice some properties of the integral transform
$\mathrm{Im}\tilde{G}_{h}(\Delta\pm\Lambda)$ (see Eq. (31)). The
characteristics of the kernel is encoded in coefficients $c_{k}$. In this way,
the Chebyshev moments $\mu_{k}$ have to be calculated only once for any kernel
to be used. This is an important feature, because their computation is much
more expensive than the post-processing stage (i.e. constructing histograms).
Moreover, the integral of Eq. (31) can be performed analytically for the
Gaussian kernel, which speeds up the calculation and does not introduce any
additional numerical errors.

###  III.3 Coupled-cluster theory

The moments of the Chebyshev expansion $\mu_{k}$ in Eq. (28) have to be
calculated in a many-body framework. In this work we employ the spherical
coupled-cluster method Hagen  _et al._ (2014), which can accurately describe
ground- and excited state properties of nuclei in the neighbourhood of closed
(sub-)shell nuclei. The method starts from a spherical Hartree-Fock reference
state $|\Psi_{\mathrm{HF}}\rangle$ and includes correlations with an
exponential ansatz

$|0\rangle=e^{T}|\Psi_{\mathrm{HF}}\rangle\,.$ (34)

Here the cluster operator $T$ is built of 1p-1h, 2p-2h, $\dots$ excitations

$T=\sum_{i,a}t^{a}_{i}a^{\dagger}_{a}a_{i}+\frac{1}{4}\sum_{ijab}t^{ab}_{ij}a_{a}^{\dagger}a_{b}^{\dagger}a_{i}a_{j}+...\,,$
(35)

and is truncated at a certain level. In this work we truncate $T$ at the 2p-2h
excitation level, which is known as the coupled-cluster singles-and-doubles
(CCSD) approximation.222In Eq. (35) index $a$, $b$ iterates over particle
states, while $i$, $j$ over hole states. The amplitudes $t$ are obtained by
solving a large set of coupled non-linear equations, which are subsequently
used in the construction of the similarity transformed Hamiltonian and
creation/annihilation operators

$\overline{H}=e^{-T}\hat{H}e^{T},\ \ \ \
\overline{a_{\alpha}}=e^{-T}a_{\alpha}e^{T},\ \ \ \
\overline{a^{\dagger}_{\alpha}}=e^{-T}a^{\dagger}_{\alpha}e^{T}\\\ .$ (36)

For our problem, the initial states are built as

$\begin{split}&|\Psi_{0}\rangle=\overline{a_{\alpha}}|\Psi_{\mathrm{HF}}\rangle\,,\
\ \ \ \
\langle\tilde{\Psi}_{0}|=\langle\Psi_{\mathrm{HF}}|\overline{a_{\beta}^{\dagger}}\,.\end{split}$
(37)

The calculation of Chebyshev moments follows Eq. (28), which requires
iterating

$\begin{split}&|\Psi_{k}\rangle=\overline{H}|\Psi_{k-1}\rangle\,,\ \ \ \ \
\langle\tilde{\Psi}_{k}|=\langle\tilde{\Psi}_{k-1}|\,\overline{H}\end{split}$
(38)

The action of the Hamiltonian can be accumulated either on the right state, or
on the left, or distributed between them. This allows for a numerical check of
the procedure.

In our calculation we use the NNLOsat nucleon-nucleon and three-nucleon
interaction, which was adjusted to the binding energy and charge radii of
light nuclei and selected oxygen isotopes Ekström  _et al._ (2015).
Furthermore, we approximate the three-nucleon interaction at the normal-
ordered two-body level which has been shown to be accurate for light- and
medium mass nuclei Hagen  _et al._ (2007); Roth  _et al._ (2012). We note that
this approximation breaks translational invariance of the Hamiltonian, and
impacts the computation of intrinsic observables in light nuclei Djärv  _et
al._ (2021). The results for various observables are converged in a model
space of 15 oscillator shells ($N_{\mathrm{max}}=14$) using the oscillator
spacing $\hbar\omega=16$ MeV. The three-nucleon interaction had an additional
energy cut on allowed configurations given by $E_{3\rm max}=16$ MeV.

##  IV Results

Before the analysis of the spectral function itself, we benchmark our result
for the momentum distribution, which is directly derived from the spectral
function:

$n(\mathbf{p})=\int
dE\,S(\mathbf{p},E)=\sum_{\alpha,\beta}\langle\mathbf{p}|\alpha\rangle\langle\mathbf{p}|\beta\rangle^{\dagger}\langle
0|a^{\dagger}_{\beta}a_{\alpha}|0\rangle\,.$ (39)

Figure 2: Comparison of laboratory and intrinsic momentum distribution in 4He
(solid and dashed black lines). The grey band corresponds to an uncertainty of
our procedure of removing the center of mass wavefunction. Results of SCGF
Rocco and Barbieri (2018) are shown for the laboratory (dashed line) and
intrinsic momentum (dashed-dotted line). All the calculations were done with
NNLOsat interaction and the same model space.

As coupled-cluster computations are performed in the laboratory system, one
has to extract the intrinsic spatial density (or intrinsic momentum
distribution) from the corresponding laboratory distributions. Within the
coupled-cluster and in-medium similarity renormalization group frameworks the
nuclear ground-state has been shown to factorize with a very good precision
into an intrinsic wavefunction and a Gaussian center-of-mass Hagen  _et al._
(2009a, 2010a); Parzuchowski  _et al._ (2017) when the kinetic energy of the
center of mass is removed from the Hamiltonian. We note that this
factorization was demonstrated in the coupled-cluster approach for a two-body
Hamiltonian, while the in-medium similarity renormalization group approach
showed that this factorization holds when applying a three-nucleon interaction
in the normal-ordered two-body approximation in nuclei as light as 14C
Parzuchowski  _et al._ (2017). Since the effect of breaking translational
invariance was found to be small on the binding energy and radius of 4He
Ekström  _et al._ (2015), we assume that factorization holds also for 4He in
the coupled-cluster theory. Assuming that the center of mass is a Gaussians,
the extraction of intrinsic momentum distribution involves a deconvolution via
Fourier transforms, and details are presented in the Appendix A. Because of
numerical reasons—varying the cut-off in the Fourier transform—the low-
momentum region is affected by a few-percent uncertainty.

Figure 2 shows the intrinsic and laboratory proton momentum distributions
computed within CCSD. The difference is clearly visible for low and
intermediate momenta up to $k\approx 0.7\mathrm{fm}^{-1}$. We compare our
results with those from the SCGF method. While the results coincide for the
laboratory momentum distribution, there are visible differences to the
intrinsic CCSD momentum density. We mostly ascribe them to two very different
strategies of the center-of-mass removal. In our case the method is
straightforward and relatively simple, while the procedure employed in Ref.
Rocco and Barbieri (2018) consists of two steps. Firstly the SCGF result is
approximated via an optimized reference state, then the center-of-mass
component is removed from the wavefunction using Monte Carlo Metropolis
sampling. We note that the three-nucleon interaction is approximated slightly
different in the SCGF approach Carbone  _et al._ (2013), and may therefore
impact the comparison with our approach for intrinsic observables in 4He. We
speculate that this difference will lead to some discrepancies in the cross-
section predictions. However, the low-momentum region plays a minor role since
the hadron tensor is weighted by $p^{2}dp$ [see Eq. (16)].

###  IV.1 Spectral function

Benchmarking of the momentum distribution $n(\mathbf{p})$ allows us to
validate the momentum dependence of the spectral function $S(\mathbf{p},E)$
and compare it with a previous calculation of Ref. Rocco and Barbieri (2018).
However, the spectral function energy dependence requires a more careful
analysis. The energy distribution, driven by $\mathrm{Im}G(\alpha,\beta,E)$ is
obtained via the integral transform expanded into Chebyshev polynomials. These
are calculated according to recursive relations of Eq. (28) iterating the
action of the Hamiltonian on the initial pivot state

$|\Psi_{n}\rangle=\hat{H}^{n}|\Psi_{0}\rangle=\hat{H}^{n}\,a_{\alpha}|0\rangle\
.$ (40)

Two remarks are in order. The initial state $|\Psi_{0}\rangle$ is composed of
$(A-1)$ nucleons and, to be consistent, the Hamiltonian applied in the
iteration should be changed accordingly. Otherwise the energy conservation of
Eq. (17), $E_{0}-E_{\Psi}$ would be shifted since $E_{\Psi}$ is the excitation
of the $(A-1)$ system. Additionally, $|\Psi_{0}\rangle$ might contain spurious
center-of-mass excitations which should be detected and removed. The first
point was already discussed in Ref. Hagen  _et al._ (2010b) and a method to
account for this inconsistency was proposed. It consists in performing the
calculation of the ground state and excitation energies both for $A$ nucleons
($E_{0}$ and $\omega=E_{0}-E_{\psi}$ accordingly) and $(A-1)$ ($E_{0}^{*}$ and
$\omega^{*}=E_{0}^{*}-E_{\psi}^{*}$ accordingly), so that we can calculate

$E_{0}-E_{\Psi}^{*}=E_{0}-E_{0}^{*}+\omega^{*}\ .$ (41)

We have checked that the difference between this value and $E_{0}-E_{\Psi}$ is
around $1.5$ MeV. For the purposes of the spectral function—which is a valid
approximation for the momentum transfer of several hundreds of MeV and energy
transfers of tens of MeV—this difference is not drastic. It will be also
partially taken into account since we consider spectral function in form of a
histogram, whose binning will be larger than $1.5$ MeV.

The disentanglement of the center of mass from the physical excitations in
$\mathrm{Im}G(\alpha,\beta,E)$ is complicated. However, the spectral function
of 4He has a simple structure. It is dominated by a single peak whose position
corresponds to the energy difference between the ground states of 4He and 3H
(in the case of the proton spectral function) or 4He and 3He (in the case of
the neutron spectral function). In Fig. 3, we show how the convergence of
$\mathrm{Im}G(\alpha,\beta,E)$ depends on $N_{\mathrm{max}}$ for the specific
case of quantum numbers $\alpha=\beta$, radial quantum number $n=0$ and
orbital angular momentum $l=0$. The dominant peak at around $23.5$ MeV is
converged already for $N_{\mathrm{max}}=14$. The smaller excitations visible
at higher energies play a negligible role. Their integrated strength does not
exceed $2\%$. We also observe some very small contribution of states with
negative strengths, which can be treated as unphysical excitations. We remove
them from the final spectral function.

The analysis for 4He shows that our treatment, although introducing some
approximations, still gives reasonable spectral functions within a few percent
of uncertainty. We were able to remove the center-of-mass contamination from
the momentum-dependent part of the spectral function, and performed various
checks to make sure that the center-of-mass excitation does not strongly
affect the energy distribution. For heavier nuclei the situation is known to
be better, since center-of-mass effects scale as $1/A$.

Figure 3: $\text{Im}G_{h}(E,\alpha,\beta)$ for $\alpha=\beta$ of 4He
corresponding to a proton with the radial quantum number $n=0$, and orbital
angular momentum $l=0$.

To obtain the final spectral function within the ChEK method we need to know
the scaling factors, see Eq. (26), which we estimate through the Lanczos
algorithm, $E_{min}=0$ and $E_{max}=200$ MeV. The reconstruction using
histograms requires also to set the bin width, and $3$ MeV is a resolution
that encompasses some of the earlier discussed uncertainties. Next, the
Gaussian kernel width is chosen such that $\Lambda$ and $\Sigma$ are kept
small (see Eq. (23)), while the number of Chebyshev moments required is still
manageable. In the results presented in this work we used the Gaussian width
$\lambda=0.12$ MeV and $\Lambda=0.5$ MeV. According to the results of Ref.
Sobczyk and Roggero (2021) the number of required moments would exceed
$N=10000$ to keep the truncation error below $1\%$. While this number seems
large it is still achievable within the CCSD approximation. Still, as was also
pointed out in Ref. Sobczyk and Roggero (2021), the current bound on the
truncation error is greatly overestimated. We have checked numerically that
results are converged already with $N=6000$ moments. With this choice of
parameters we obtain a histogram spectral function whose errors are
negligible.

###  IV.2 Electron-nucleus scattering

Figure 4: Electron scattering on 4He for different kinematics which correspond
to the momentum transfers $|\mathbf{q}|\approx 270-670$ MeV. We show results
obtained with CCSD before removing the center-of-mass contamination (dashed
line) and afterwards, as explained in the main text (dashed-dotted line).
Experimental data was taken from Ref. O’Connell  _et al._ (1987); Zghiche  _et
al._ (1994); Sealock  _et al._ (1989); Day  _et al._ (1993).

We now turn our attention to the electron scattering off 4He. In Fig. 4, we
show results for the cross section in various kinematics, for the spectral
function before and after removal of the center-of-mass contamination. The
final results, “CCSD intrinsic”, has been obtained using the intrinsic
momentum distribution. They predict more strength at the quasi-elasic peak
with respect to the “CCSD lab” result. This trend is consistent with the
findings of Ref. Rocco and Barbieri (2018). The uncertainty of $n(\mathbf{p})$
at low momenta as well as the negligible reconstruction errors coming from the
ChEK method do not affect the cross-section results. It is also interesting to
notice that the impulse approximation indeed works better with increasing
momentum transfer $|\mathbf{q}|$. For the values $|\mathbf{q}|\approx 300-400$
MeV the spectral function overestimates the data and predicts a shifted quasi-
elastic peak.

A direct comparison of our results with Ref. Rocco and Barbieri (2018) shows
an overall good agreement. While the results before the center-of-mass removal
are almost identical, the predicted cross section using the intrinsic spectral
function is slightly different, as can be seen in Fig 5. We have chosen this
low momentum transfer kinematics, because the nuclear effects might play a
more important role and differences between the CCSD and SCGF should be more
pronounced. There are several sources of discrepancies. First, in the
conservation of energy of Eq. (II.1) we take into account the kinetic energy
of the recoiled nucleus, which for 4He amounts to 7-9 MeV for the Fermi
momentum. Second, we use a different approach to remove the center-of-mass.
Third, the spectral functions are obtained using two different many-body
methods and approximations therein.

Figure 5: Comparison of CCSD results with SCFG Rocco and Barbieri (2018) on
4He after removing the center-of-mass contamination.

##  V Conclusion and outlook

We have presented an ab initio calculation of the spectral functions for 4He
based on the coupled cluster theory combined with the ChEK method for the
reconstruction of the spectral properties of a many-body system. Within this
approach, and for a given resolution, we were able to assess the uncertainty
of our calculation. For 4He we obtained an almost negligible error, however we
expect it to be larger when we move to medium-mass nuclei. This work paves the
way for further explorations of the ChEK method in nuclear physics. On the one
hand, one can study other nuclear responses, especially where the standard
inversion procedures are unable to give stable results. On the other hand, the
method provides error bounds, because it does not require any ansatz about the
properties of the response (e.g. the threshold energy needed in the LIT
inversion). This way we can achieve a stringent control over the uncertainty
bound of an observable of interest.

We compared our predictions for the electron-nucleus scattering in the quasi-
elastic regime with available data and found agreement. We were able to scan a
large range of momentum transfers, and observed that the impulse approximation
improves with growing momentum. Our results point to possible directions of
further investigation. In view of the planned DUNE and T2HK experiments, which
will benefit from reliable cross-section models, another important step is the
calculation of spectral function for 16O and 40Ar. Both of these nuclei is
within the reach of the coupled-cluster approach. Furthermore, for low and
intermediate momentum transfers the impulse approximation picture becomes less
reliable and the full inclusion of final state interactions is therefore
desirable. This transition region requires more theoretical studies as well as
experimental data. Recent developments within coupled-cluster theory allow us
to lead a consistent analysis based on the same many-body method, nuclear
dynamics and truncations, employing both LIT-CC and spectral functions.
Lastly, our approach does not presently account for two-body currents. Their
role will also be a topic of our future work.

###### Acknowledgements.

J.E.S. acknowledges the support of the Humboldt Foundation through a Humboldt
Research Fellowship for Postdoctoral Researchers. This work was supported in
part by the Deutsche Forschungsgemeinschaft (DFG) through the Cluster of
Excellence “Precision Physics, Fundamental Interactions, and Structure of
Matter” (PRISMA+ EXC 2118/1) funded by the DFG within the German Excellence
Strategy (Project ID 39083149). This material is based upon work supported by
the U.S. Department of Energy, Office of Science, Office of Nuclear Physics
under award numbers DE-FG02-96ER40963 and DE-SC0018223 (NUCLEI SciDAC-4
collaboration), and contract No. DE-AC05-00OR22725 with UT-Battelle, LLC (Oak
Ridge National Laboratory). Computer time was provided by the Innovative and
Novel Computational Impact on Theory and Experiment (INCITE) program and by
the supercomputer Mogon at Johannes Gutenberg Universität Mainz. This research
used resources of the Oak Ridge Leadership Computing Facility located at Oak
Ridge National Laboratory, which is supported by the Office of Science of the
Department of Energy under contract No. DE-AC05-00OR22725.

##  Appendix A Removal of the center of mass

We follow Ref. Giraud (2008) and consider an $A$-body system with coordinates
$\mathbf{r}_{1},\ldots,\mathbf{r}_{A}$ and corresponding momenta
$\mathbf{p}_{1},\ldots,\mathbf{p}_{A}$ in the laboratory system. The center-
of-mass and relative coordinates are

$\displaystyle\mathbf{R}$ $\displaystyle=$ $\displaystyle{1\over
A}\left(\mathbf{r}_{1}+\ldots+\mathbf{r}_{A}\right)\ ,$
$\displaystyle\bm{\xi}_{1}$ $\displaystyle=$
$\displaystyle\mathbf{r}_{2}-\mathbf{r}_{1}\ ,$ $\displaystyle\bm{\xi}_{2}$
$\displaystyle=$
$\displaystyle\mathbf{r}_{3}-\frac{\mathbf{r}_{1}+\mathbf{r}_{2}}{2}\ ,$
$\displaystyle\vdots$ $\displaystyle\bm{\xi}_{n}$ $\displaystyle=$
$\displaystyle\mathbf{r}_{n+1}-\frac{\mathbf{r}_{1}+\mathbf{r}_{2}+\ldots+\mathbf{r}_{n}}{n}\
,$ $\displaystyle\vdots$ $\displaystyle\bm{\xi}_{A-1}$ $\displaystyle=$
$\displaystyle\mathbf{r}_{A}-\frac{\mathbf{r}_{1}+\mathbf{r}_{2}+\ldots+\mathbf{r}_{A-1}}{A-1}\
,$ (42)

and the corresponding canonical momenta are

$\displaystyle\mathbf{P}$ $\displaystyle=$
$\displaystyle\mathbf{p}_{1}+\ldots+\mathbf{p}_{A}\ ,$
$\displaystyle\bm{\pi}_{1}$ $\displaystyle=$
$\displaystyle\frac{\mathbf{p}_{2}-\mathbf{r}_{1}}{2}\ ,$
$\displaystyle\bm{\pi}_{2}$ $\displaystyle=$ $\displaystyle{2\over
3}\mathbf{p}_{3}-\frac{\mathbf{p}_{1}+\mathbf{p}_{2}}{3}\ ,$
$\displaystyle\vdots$ $\displaystyle\bm{\pi}_{n}$ $\displaystyle=$
$\displaystyle{n\over
n+1}\mathbf{p}_{n+1}-\frac{\mathbf{p}_{1}+\mathbf{p}_{2}+\ldots+\mathbf{p}_{n}}{n+1}\
,$ $\displaystyle\vdots$ $\displaystyle\bm{\pi}_{A-1}$ $\displaystyle=$
$\displaystyle{A-1\over
A}\mathbf{p}_{A}-\frac{\mathbf{p}_{1}+\mathbf{p}_{2}+\ldots+\mathbf{p}_{A-1}}{A}\
.$ (43)

Two comments are in order. First, the transformation between laboratory and
center-of-mass coordinates has a unit Jacobian. Second, we see that
$\bm{\xi}_{n}$ is the position of particle $(n+1)$ with respect to the center
of mass of the previous $n$ particles, while $\bm{\pi}_{n}$ is the momentum of
particle $(n+1)$ relative to the average momentum of the $(n+1)$ particle
system. We note that $\mathbf{r}_{i}-\mathbf{R}$ and
$\mathbf{p}_{i}-\mathbf{P}/A$ are intrinsic positions and momenta,
respectively, for any $i=1,\ldots,A$. A particular convenient choice is $i=A$,
because

$\displaystyle\mathbf{r}_{A}-\mathbf{R}$ $\displaystyle=$
$\displaystyle{A-1\over A}\bm{\xi}_{A-1}\ ,$ (44)
$\displaystyle\mathbf{p}_{A}-{\mathbf{P}\over A}$ $\displaystyle=$
$\displaystyle\bm{\pi}_{A-1}$ (45)

can be expressed in terms of a single relative position and momentum,
respectively.

Let us now assume that the ground state
$|\Psi\rangle=|\Phi\rangle|\psi\rangle$ factorizes into a center-of-mass state
$|\Phi\rangle$ and an intrinsic state $|\psi\rangle$. The laboratory density
in position space is

$\displaystyle\rho(\mathbf{x})$ $\displaystyle=$ $\displaystyle A\int{\rm
d}^{3}\mathbf{r}_{1}\cdots{\rm
d}^{3}\mathbf{r}_{A}\delta(\mathbf{x}-\mathbf{r}_{A})|\Psi(\mathbf{r}_{1},\ldots,\mathbf{r}_{A})|^{2}$
(46) $\displaystyle=$ $\displaystyle A\int{\rm d}^{3}\mathbf{r}_{1}\cdots{\rm
d}^{3}\mathbf{r}_{A-1}|\Psi(\mathbf{r}_{1},\ldots,\mathbf{r}_{A-1},\mathbf{x})|^{2}\,.$

Based on Eq. (44), the intrinsic one-body density is

$\displaystyle\sigma(\mathbf{x})$ $\displaystyle=$ $\displaystyle A\int{\rm
d}^{3}\bm{\xi}_{1}\cdots{\rm
d}^{3}\bm{\xi}_{A-1}\delta\left(\mathbf{x}-{A-1\over A}\bm{\xi}_{A-1}\right)$
$\displaystyle\times|\psi(\bm{\xi}_{1},\ldots,\bm{\xi}_{A-1})|^{2}$
$\displaystyle=$ $\displaystyle A\left({A\over A-1}\right)^{3}$
$\displaystyle\times\int{\rm d}^{3}\bm{\xi}_{1}\cdots{\rm
d}^{3}\bm{\xi}_{A-2}\left|\psi\left(\bm{\xi}_{1},\ldots,\bm{\xi}_{A-2},{A\over
A-1}\mathbf{x}\right)\right|^{2}\,.$

To establish the relation between the intrinsic and laboratory densities we
use
$\delta(\mathbf{x}-\mathbf{r}_{A})=\delta(\mathbf{x}-\mathbf{R}-\bm{\xi}_{A-1}(A-1)/A)$
and rewrite Eq. (46) as

$\displaystyle\rho(\mathbf{x})$ $\displaystyle=$ $\displaystyle A\int{\rm
d}^{3}\bm{\xi}_{1}\cdots{\rm d}^{3}\bm{\xi}_{A-1}{\rm
d}^{3}\mathbf{R}\delta\left(\mathbf{x}-\mathbf{R}-{A-1\over
A}\bm{\xi}_{A-1}\right)$ (48)
$\displaystyle\times|\Phi(\mathbf{R})|^{2}|\psi(\bm{\xi}_{1},\ldots,\bm{\xi}_{A-1})|^{2}$
$\displaystyle=$ $\displaystyle\int{\rm
d}^{3}\mathbf{R}|\Phi(\mathbf{R})|^{2}\sigma(\mathbf{x}-\mathbf{R})\ .$

Thus, the laboratory density is a convolution of the center-of-mass density
and the intrinsic density. In coupled-cluster and IMSRG calculations, the
center-of-mass wave function is a Gaussian (to a very good approximation)
Hagen  _et al._ (2009b, 2010a); Jansen (2013); Morris  _et al._ (2015);
Parzuchowski  _et al._ (2017), i.e.

$|\Phi(\mathbf{R})|^{2}=\pi^{-3/2}b^{-3}e^{-{\mathbf{R}^{2}\over b^{2}}}\ ,$
(49)

and its Fourier transform is

$(2\pi)^{-3/2}e^{-{b^{2}\over 4}\mathbf{P}^{2}}\ .$ (50)

Thus, the intrinsic density can be obtained by dividing the Fourier transforms
of the laboratory and center-of-mass wave functions and performing the inverse
Fourier transform of that quotient.

Let us now turn to momentum space. The laboratory momentum density is

$\displaystyle\rho(\mathbf{p})$ $\displaystyle=$ $\displaystyle A\int{\rm
d}^{3}\mathbf{p}_{1}\cdots{\rm
d}^{3}\mathbf{p}_{A}\delta(\mathbf{p}-\mathbf{p}_{A})|\Psi(\mathbf{p}_{1},\ldots,\mathbf{p}_{A})|^{2}$
(51) $\displaystyle=$ $\displaystyle A\int{\rm d}^{3}\mathbf{p}_{1}\cdots{\rm
d}^{3}\mathbf{p}_{A-1}|\Psi(\mathbf{p}_{1},\ldots,\mathbf{p}_{A-1},\mathbf{p})|^{2}\,.$

Based on Eq. (45), the intrinsic one-body momentum density is

$\displaystyle\sigma(\mathbf{p})$ $\displaystyle=$ $\displaystyle A\int{\rm
d}^{3}\bm{\pi}_{1}\cdots{\rm
d}^{3}\bm{\pi}_{A-1}\delta\left(\mathbf{p}-\bm{\xi}_{A-1}\right)$
$\displaystyle\times|\psi(\bm{\pi}_{1},\ldots,\bm{\pi}_{A-1})|^{2}$
$\displaystyle=$ $\displaystyle A\int{\rm d}^{3}\bm{\pi}_{1}\cdots{\rm
d}^{3}\bm{\pi}_{A-2}\left|\psi\left(\bm{\pi}_{1},\ldots,\bm{\pi}_{A-2},\mathbf{p}\right)\right|^{2}\
.$

To establish the relation between the intrinsic and laboratory momentum
densities we use
$\delta(\mathbf{p}-\mathbf{p}_{A})=\delta(\mathbf{p}-\bm{\pi}_{A-1}-\mathbf{P}/A)$
and rewrite Eq. (51) as

$\displaystyle\rho(\mathbf{p})$ $\displaystyle=$ $\displaystyle A\int{\rm
d}^{3}\bm{\pi}_{1}\cdots{\rm d}^{3}\bm{\pi}_{A-1}{\rm
d}^{3}\mathbf{P}\delta\left(\mathbf{p}-\bm{\pi}_{A-1}-{\mathbf{P}\over
A}\right)$ (53)
$\displaystyle\times|\Phi(\mathbf{P})|^{2}|\psi(\bm{\pi}_{1},\ldots,\bm{\pi}_{A-1})|^{2}$
$\displaystyle=$ $\displaystyle\int{\rm
d}^{3}\mathbf{P}|\Phi(\mathbf{P})|^{2}\sigma\left(\mathbf{p}-{\mathbf{P}\over
A}\right)\ .$

The substitition $\mathbf{K}=\mathbf{P}/A$ then yields

$\displaystyle\rho(\mathbf{p})=A^{3}\int{\rm
d}^{3}\mathbf{K}|\Phi(A\mathbf{K})|^{2}\sigma\left(\mathbf{p}-\mathbf{K}\right)\
.$ (54)

Thus, the laboratory density is a convolution of the center-of-mass density
(at $A$ times its argument) and the intrinsic density. Again, the center-of-
mass wave function is to a good approximation a Gaussian in momentum space,
and the deconvolution can be performed. We have

$A^{3}|\Phi(A\mathbf{P})|^{2}=\pi^{-3/2}b^{3}A^{3}e^{-b^{2}A^{2}\mathbf{P}^{2}}\
,$ (55)

and its Fourier transform is

$(2\pi)^{-3/2}e^{-{1\over 4}{\mathbf{R}^{2}\over A^{2}b^{2}}}\ .$ (56)

In the Gaussians, we employed the oscillator length

$b={\sqrt{\hbar\over M\tilde{\omega}}}\ .$ (57)

Here, $M=Am$ is the total mass in terms of the nucleon mass $m$, while
$\tilde{\omega}$ is the frequency of the Gaussian; this parameter is
independent of the oscillator basis, and we have $\hbar\tilde{\omega}\approx
24$ MeV for a light nucleus such as 4He and $\hbar\tilde{\omega}\approx 16-20$
MeV for 16OHagen  _et al._ (2009b, 2010a). Thus $b\approx 0.7A^{-1/2}$ fm, and
the Gaussian (55) approximates a $\delta$-function for $A\gg 1$. We therefore
expect that the difference between the laboratory and intrinsic density in
momentum space becomes exponentially small for $A^{1/2}\gg 1$.

## References

  * Bacca and Pastore (2014) Sonia Bacca and Saori Pastore, “Electromagnetic reactions on light nuclei,” Journal of Physics G: Nuclear and Particle Physics 41, 123002 (2014).
  * Alvarez-Ruso  _et al._ (2018) L. Alvarez-Ruso _et al._ (NuSTEC), “NuSTEC White Paper: Status and challenges of neutrino–nucleus scattering,” Prog. Part. Nucl. Phys. 100, 1–68 (2018), arXiv:1706.03621 [hep-ph] .
  * Ankowski  _et al._ (2022) A. M. Ankowski _et al._ , “Electron Scattering and Neutrino Physics,” (2022), arXiv:2203.06853 [hep-ex] .
  * Ruso  _et al._ (2022) L. Alvarez Ruso _et al._ , “Theoretical tools for neutrino scattering: interplay between lattice QCD, EFTs, nuclear physics, phenomenology, and neutrino event generators,” (2022), arXiv:2203.09030 [hep-ph] .
  * Bacca  _et al._ (2013) S. Bacca, N. Barnea, G. Hagen, G. Orlandini, and T. Papenbrock, “First principles description of the giant dipole resonance in ${}^{16}\mathbf{O}$,” Phys. Rev. Lett. 111, 122502 (2013).
  * Acharya and Bacca (2020) Bijaya Acharya and Sonia Bacca, “Neutrino-deuteron scattering: Uncertainty quantification and new $L_{1,A}$ constraints,” Phys. Rev. C 101, 015505 (2020), arXiv:1911.12659 [nucl-th] .
  * Sobczyk  _et al._ (2020) J. E. Sobczyk, B. Acharya, S. Bacca, and G. Hagen, “Coulomb sum rule for 4He and 16O from coupled-cluster theory,” Phys. Rev. C 102, 064312 (2020), arXiv:2009.01761 [nucl-th] .
  * Sobczyk  _et al._ (2021) J. E. Sobczyk, B. Acharya, S. Bacca, and G. Hagen, “Ab initio computation of the longitudinal response function in 40Ca,” Phys. Rev. Lett. 127, 072501 (2021), arXiv:2103.06786 [nucl-th] .
  * Nieves  _et al._ (2011) J. Nieves, I. Ruiz Simo, and M. J. Vicente Vacas, “Inclusive Charged–Current Neutrino–Nucleus Reactions,” Phys. Rev. C 83, 045501 (2011), arXiv:1102.2777 [hep-ph] .
  * Martinez  _et al._ (2006) M. C. Martinez, P. Lava, N. Jachowicz, J. Ryckebusch, K. Vantournhout, and J. M. Udias, “Relativistic models for quasi-elastic neutrino scattering,” Phys. Rev. C 73, 024607 (2006), arXiv:nucl-th/0505008 .
  * Amaro  _et al._ (2005) Jose Enrique Amaro, M. B. Barbaro, J. A. Caballero, T. W. Donnelly, A. Molinari, and I. Sick, “Using electron scattering superscaling to predict charge-changing neutrino cross sections in nuclei,” Phys. Rev. C 71, 015501 (2005), arXiv:nucl-th/0409078 .
  * Leitner  _et al._ (2009) T. Leitner, O. Buss, L. Alvarez-Ruso, and U. Mosel, “Electron- and neutrino-nucleus scattering from the quasielastic to the resonance region,” Phys. Rev. C 79, 034601 (2009), arXiv:0812.0587 [nucl-th] .
  * Benhar  _et al._ (2005) Omar Benhar, Nicola Farina, Hiroki Nakamura, Makoto Sakuda, and Ryoichi Seki, “Electron- and neutrino-nucleus scattering in the impulse approximation regime,” Phys. Rev. D 72, 053005 (2005), arXiv:hep-ph/0506116 .
  * Benhar  _et al._ (1994) O. Benhar, A. Fabrocini, S. Fantoni, and I. Sick, “Spectral function of finite nuclei and scattering of GeV electrons,” Nucl. Phys. A 579, 493–517 (1994).
  * Nieves and Sobczyk (2017) Juan Nieves and Joanna Ewa Sobczyk, “In medium dispersion relation effects in nuclear inclusive reactions at intermediate and low energies,” Annals Phys. 383, 455–496 (2017), arXiv:1701.03628 [nucl-th] .
  * Buss  _et al._ (2007) O. Buss, T. Leitner, U. Mosel, and L. Alvarez-Ruso, “The Influence of the nuclear medium on inclusive electron and neutrino scattering off nuclei,” Phys. Rev. C 76, 035502 (2007), arXiv:0707.0232 [nucl-th] .
  * Rocco and Barbieri (2018) N. Rocco and C. Barbieri, “Inclusive electron-nucleus cross section within the Self Consistent Green’s Function approach,” Phys. Rev. C 98, 025501 (2018), arXiv:1803.00825 [nucl-th] .
  * Barbieri  _et al._ (2019) C. Barbieri, N. Rocco, and V. Somà, “Lepton Scattering from 40Ar and Ti in the Quasielastic Peak Region,” Phys. Rev. C 100, 062501 (2019), arXiv:1907.01122 [nucl-th] .
  * Efros  _et al._ (1998) Victor D. Efros, Winfried Leidemann, and Giuseppina Orlandini, “Exact He-4 spectral function in a semirealistic N N potential model,” Phys. Rev. C 58, 582–585 (1998), arXiv:nucl-th/9805002 .
  * Efros  _et al._ (1994) Victor D. Efros, Winfried Leidemann, and Giuseppina Orlandini, “Response functions from integral transforms with a Lorentz kernel,” Phys. Lett. B 338, 130–133 (1994), arXiv:nucl-th/9409004 .
  * Efros  _et al._ (2007) V D Efros, W Leidemann, G Orlandini, and N Barnea, “The lorentz integral transform (lit) method and its applications to perturbation-induced reactions,” Journal of Physics G: Nuclear and Particle Physics 34, R459 (2007).
  * Carlson and Schiavilla (1992) J. Carlson and R. Schiavilla, “Euclidean proton response in light nuclei,” Phys. Rev. Lett. 68, 3682–3685 (1992).
  * Carlson  _et al._ (2015) J. Carlson, S. Gandolfi, F. Pederiva, Steven C. Pieper, R. Schiavilla, K. E. Schmidt, and R. B. Wiringa, “Quantum monte carlo methods for nuclear physics,” Rev. Mod. Phys. 87, 1067–1118 (2015).
  * Lovato  _et al._ (2016) A. Lovato, S. Gandolfi, J. Carlson, Steven C. Pieper, and R. Schiavilla, “Electromagnetic response of ${}^{12}\mathrm{C}$: A first-principles calculation,” Phys. Rev. Lett. 117, 082501 (2016).
  * Bacca  _et al._ (2014) S. Bacca, N. Barnea, G. Hagen, M. Miorelli, G. Orlandini, and T. Papenbrock, “Giant and pigmy dipole resonances in ${}^{4}\mathrm{He}$, ${}^{16,22}\mathrm{O}$, and ${}^{40}\mathrm{Ca}$ from chiral nucleon-nucleon interactions,” Phys. Rev. C 90, 064619 (2014).
  * Lovato  _et al._ (2020) A. Lovato, J. Carlson, S. Gandolfi, N. Rocco, and R. Schiavilla, “Ab initio study of $({\nu}_{\ell},{\ell}^{-})$ and $({\overline{\nu}}_{\ell},{\ell}^{+})$ inclusive scattering in ${}^{12}\mathrm{C}$: Confronting the miniboone and t2k ccqe data,” Phys. Rev. X 10, 031068 (2020).
  * Bacca  _et al._ (2002) Sonia Bacca, Mario Andrea Marchisio, Nir Barnea, Winfried Leidemann, and Giuseppina Orlandini, “Microscopic calculation of six-body inelastic reactions with complete final state interaction: Photoabsorption of 6he and 6li,” Phys. Rev. Lett. 89, 052502 (2002).
  * Bacca  _et al._ (2004) Sonia Bacca, Nir Barnea, Winfried Leidemann, and Giuseppina Orlandini, “Effect of $p$-wave interaction in ${}^{6}\mathrm{He}$ and ${}^{6}\mathrm{Li}$ photoabsorption,” Phys. Rev. C 69, 057001 (2004).
  * Raghavan  _et al._ (2021) Krishnan Raghavan, Prasanna Balaprakash, Alessandro Lovato, Noemi Rocco, and Stefan M. Wild, “Machine-learning-based inversion of nuclear responses,” Phys. Rev. C 103, 035502 (2021).
  * Roggero (2020) Alessandro Roggero, “Spectral density estimation with the Gaussian Integral Transform,” Phys. Rev. A 102, 022409 (2020), arXiv:2004.04889 [quant-ph] .
  * Sobczyk and Roggero (2021) Joanna E. Sobczyk and Alessandro Roggero, “Spectral density reconstruction with Chebyshev polynomials,” (2021), arXiv:2110.02108 [nucl-th] .
  * Lanczos (1950) Cornelius Lanczos, “An iteration method for the solution of the eigenvalue problem of linear differential and integral operators,” J. Res. Natl. Bur. Stand. B 45, 255–282 (1950).
  * Weiße  _et al._ (2006) Alexander Weiße, Gerhard Wellein, Andreas Alvermann, and Holger Fehske, “The kernel polynomial method,” Reviews of Modern Physics 78, 275–306 (2006).
  * Bradford  _et al._ (2006) R. Bradford, A. Bodek, Howard Scott Budd, and J. Arrington, “A New parameterization of the nucleon elastic form-factors,” Nucl. Phys. B Proc. Suppl. 159, 127–132 (2006), arXiv:hep-ex/0602017 .
  * Dickhoff and Van Neck (2005) Willem H Dickhoff and Dimitri Van Neck,  _Many-Body Theory Exposed!_ (World Scientific, 2005).
  * Tropiano  _et al._ (2021) A. J. Tropiano, S. K. Bogner, and R. J. Furnstahl, “Short-range correlation physics at low renormalization group resolution,” Phys. Rev. C 104, 034311 (2021).
  * Suzuki  _et al._ (2005) Ryusuke Suzuki, Takayuki Myo, and Kiyoshi Kato, “Level density in complex scaling method,” AIP Conf. Proc. 768, 455 (2005), arXiv:nucl-th/0502012 .
  * Carbonell  _et al._ (2014) J. Carbonell, A. Deltuva, A. C. Fonseca, and R. Lazauskas, “Bound state techniques to solve the multiparticle scattering problem,” Prog. Part. Nucl. Phys. 74, 55–80 (2014), arXiv:1310.6631 [nucl-th] .
  * Rotureau  _et al._ (2017) J. Rotureau, P. Danielewicz, G. Hagen, F. Nunes, and T. Papenbrock, “Optical potential from first principles,” Phys. Rev. C 95, 024315 (2017), arXiv:1611.04554 [nucl-th] .
  * Rotureau  _et al._ (2018) J. Rotureau, P. Danielewicz, G. Hagen, G. R. Jansen, and F. M. Nunes, “Microscopic optical potentials for calcium isotopes,” Phys. Rev. C 98, 044625 (2018), arXiv:1808.04535 [nucl-th] .
  * Hagen  _et al._ (2014) G. Hagen, T. Papenbrock, M. Hjorth-Jensen, and D. J. Dean, “Coupled-cluster computations of atomic nuclei,” Rep. Prog. Phys. 77, 096302 (2014).
  * Ekström  _et al._ (2015) A. Ekström, G. R. Jansen, K. A. Wendt, G. Hagen, T. Papenbrock, B. D. Carlsson, C. Forssén, M. Hjorth-Jensen, P. Navrátil, and W. Nazarewicz, “Accurate nuclear radii and binding energies from a chiral interaction,” Phys. Rev. C 91, 051301 (2015), arXiv:1502.04682 [nucl-th] .
  * Hagen  _et al._ (2007) G. Hagen, T. Papenbrock, D. J. Dean, A. Schwenk, A. Nogga, M. Włoch, and P. Piecuch, “Coupled-cluster theory for three-body Hamiltonians,” Phys. Rev. C 76, 034302 (2007).
  * Roth  _et al._ (2012) Robert Roth, Sven Binder, Klaus Vobig, Angelo Calci, Joachim Langhammer, and Petr Navrátil, “Medium-Mass Nuclei with Normal-Ordered Chiral $NN\mathbf{+}3N$ Interactions,” Phys. Rev. Lett. 109, 052501 (2012).
  * Djärv  _et al._ (2021) T. Djärv, A. Ekström, C. Forssén, and G. R. Jansen, “Normal-ordering approximations and translational (non)invariance,” Phys. Rev. C 104, 024324 (2021).
  * Hagen  _et al._ (2009a) G. Hagen, T. Papenbrock, and D. J. Dean, “Solution of the center-of-mass problem in nuclear structure calculations,” Phys. Rev. Lett. 103, 062503 (2009a), arXiv:0905.3167 [nucl-th] .
  * Hagen  _et al._ (2010a) G. Hagen, T. Papenbrock, D. J. Dean, and M. Hjorth-Jensen, “Ab initio coupled-cluster approach to nuclear structure with modern nucleon-nucleon interactions,” Phys. Rev. C 82, 034330 (2010a).
  * Parzuchowski  _et al._ (2017) N. M. Parzuchowski, S. R. Stroberg, P. Navrátil, H. Hergert, and S. K. Bogner, “Ab initio electromagnetic observables with the in-medium similarity renormalization group,” Phys. Rev. C 96, 034324 (2017).
  * Carbone  _et al._ (2013) Arianna Carbone, Andrea Cipollone, Carlo Barbieri, Arnau Rios, and Artur Polls, “Self-consistent green’s functions formalism with three-body interactions,” Phys. Rev. C 88, 054326 (2013).
  * Hagen  _et al._ (2010b) G. Hagen, T. Papenbrock, and M. Hjorth-Jensen, “Ab Initio Computation of the F-17 Proton Halo State and Resonances in A=17 Nuclei,” Phys. Rev. Lett. 104, 182501 (2010b), arXiv:1003.1995 [nucl-th] .
  * O’Connell  _et al._ (1987) J. S. O’Connell _et al._ , “Electromagnetic excitation of the delta resonance in nuclei,” Phys. Rev. C35, 1063 (1987).
  * Zghiche  _et al._ (1994) A. Zghiche _et al._ , “Longitudinal and transverse responses in quasielastic electron scattering from pb-208 and he-4,” Nucl. Phys. A572, 513–559 (1994).
  * Sealock  _et al._ (1989) R. M. Sealock _et al._ , “Electroexcitation of the delta (1232) in nuclei,” Phys. Rev. Lett. 62, 1350–1353 (1989).
  * Day  _et al._ (1993) D. B. Day _et al._ , “Inclusive electron nucleus scattering at high momentum transfer,” Phys. Rev. C48, 1849–1863 (1993).
  * Giraud (2008) B. G. Giraud, “Density functionals in the laboratory frame,” Phys. Rev. C 77, 014311 (2008).
  * Hagen  _et al._ (2009b) G. Hagen, T. Papenbrock, and D. J. Dean, “Solution of the center-of-mass problem in nuclear structure calculations,” Phys. Rev. Lett. 103, 062503 (2009b).
  * Jansen (2013) G. R. Jansen, “Spherical coupled-cluster theory for open-shell nuclei,” Phys. Rev. C 88, 024305 (2013).
  * Morris  _et al._ (2015) T. D. Morris, N. M. Parzuchowski, and S. K. Bogner, “Magnus expansion and in-medium similarity renormalization group,” Phys. Rev. C 92, 034331 (2015).


# Predicting Winning Captions for Weekly New Yorker Comics

Stanley Cao
Stanford University
Department of Computer Science
<EMAIL_ADDRESS>Sonny Young
Stanford University
Department of Chemistry
<EMAIL_ADDRESS>
###### Abstract

Image captioning using Vision Transformers (ViTs) represents a pivotal
convergence of computer vision and natural language processing, offering the
potential to enhance user experiences, improve accessibility, and provide
textual representations of visual data. This paper explores the application of
image captioning techniques to New Yorker cartoons, aiming to generate
captions that emulate the wit and humor of winning entries in the New Yorker
Cartoon Caption Contest. This task necessitates sophisticated visual and
linguistic processing, along with an understanding of cultural nuances and
humor. We propose several new baselines for using vision transformer encoder-
decoder models to generate captions for the New Yorker cartoon caption
contest.

##  1 Introduction

Image captioning using Vision Transformers (ViTs) is a pivotal area of
research, merging the domains of computer vision and natural language
processing. The ability to automatically generate accurate and contextually
relevant descriptions for images comes with numerous benefits, including
enhancing user experience in various digital platforms, increasing
accessibility for visually impaired users, and providing textual
representations of visual data.

In this paper, we explore the complexities of applying image captioning
techniques to New Yorker cartoons, particularly focusing on generating
captions that closely emulate the wit and humor of winning entries in the New
Yorker Cartoon Caption Contest. This task requires not only a sophisticated
combination of visual and linguistic processing, but also the ability to
caption images in a completely orthogonal way: the winning caption rarely
describes the cartoon, but rather presupposes a deep understanding of the
scene along with the concomitant cultural references. By aiming to create
models that can effectively mimic human-like humor in captioning, we
demonstrate potential advancements in the development of AI systems that have
more holistic multi-modal understanding.

Image Description |  Uncanny Description |  Entities |  Caption
---|---|---|---
Three dogs are in a field in front of some rolling hills. One of the dogs is holding a stick and winding up, as if it’s about to hit a baseball, while the other two dogs look on. |  It’s unusual that the dog is playing a sport like baseball because dogs don’t have the cognitive capacity to play that human game. |  Fetch, Dog, Anthropomorphism. |  He’s his own best friend.

Figure 1: Detailed Information including selected metadata entries for New
Yorker Cartoon Caption Contest #102

##  2 Related Work

Image captioning has so far been applied to numerous generative and
classification problems dealing with satellite data [22], pathological
diagnosis [23, 17], or real-life object detection [1, 28, 10] that primarily
deal with natural images. Much success in these tasks has been attributed to
attention models popularized by transformers [24] which embed sequential data
in a parallel and expressive manner, and more recently, ViTs which are
considered the state of the art architecture for image recognition [4].

The challenge is especially apparent when dealing with artificially
constructed images, such as cartoons or comics. These images often contain
abstract, stylized representations that may not correspond directly to real-
world objects and may embed complex narratives and humor, which are
particularly difficult for AI to parse and understand. Furthermore, ViTs are
known to require significant pre-training data, close to 100 million images
before its performance is comparable to a CNN [4]. Some have attempted to
circumvent the large dataset issue by training ViTs without using natural data
[14]. One particularly noteworthy model is LLaVA, a vision language model
(VLM) that is capable of incorporating multi-modal information in its
generated output [13]. [11] showcases that with simple modifications to the
multi-modal projection layer and further finetuning on VQA datasets, the model
is able to achieve improved performance on various downstream tasks. However,
evaluating the model beyond traditional question-answering tasks remains
largely unexplored.

Natural Language Inference (NLI) is a key area that can enhance a model’s
ability to interpret language, which is crucial for improving how AI systems
grasp concepts like humor [2], common sense [25], and overall sentence
interpretation [29]. These developments suggest a path forward for refining AI
approaches to better understand and engage with human communication nuances.

We focus on the New Yorker Cartoon Caption Contest, which requires humorous
captioning that transcends pure information content. Each week, participants
are invited to submit captions for a cartoon, often characterized by its
humorous, satirical, or ironic undertone, making the task of generating a
fitting caption not just a test of understanding visual elements, but also of
grasping subtle cultural and contextual nuances. Successfully modeling the New
Yorker Cartoon captions would demonstrate a significant leap towards AI
systems that can understand and generate human-like text responses.

Notably, [5] has explored this area of research, not only creating a dataset
of New Yorker cartoons from prior contests, but proposing three different
tasks to benchmark the capabilities revolving around multi-modal humor
understanding: matching, quality ranking, and explanation. The authors report
GPT-4 having superior performance in matching and quality ranking, achieving
84.5% accuracy in matching a true caption among a pool of 5 fake captions.
GPT-4 significantly outperformed all other methods, but it was unable to
explain jokes as well as humans. The tasks proposed in [5] are largely
classification tasks, with the exception of their proposed explanation task.
Since LLMs are unable to achieve human-level performance, these results
suggest that there is room for further exploration in the generation of
humorous captions.

##  3 Dataset

We utilized the New Yorker Caption Contest dataset available on Hugging Face
[7, 18, 21]. This dataset consists of cartoons that are used for three
different tasks to assess the multi-modal capabilities of ViTs and encoder-
decoder models. Since the task is to generate captions from the cartoon alone,
much like what is expected of humans for the caption contest, we make use of
the dataset that corresponds to the explanation task proposed in [5], which
not only contains an cartoon and winning caption pair, but also important
metadata such as human-annotated descriptions of the cartoon, an explanation
for why the caption is humorous, important entities that are relevant to the
cartoon, and leading questions that appear to motivate the winning caption.
Figure 1 shows a single data entry with selected metadata.

This dataset comprises approximately 2.6K data entries. We use the same train-
validation-test splits as provided in [7], with about 2.34k entries in the
train set and about 130 entries in the validation and test set. We do not
preprocess or normalize the cartoons in advance, as the ViTs used in this
paper contain their own custom image transformations / processing techniques
(see Section 4 for more details). As all cartoons are of different sizes, we
have varying image resolutions, though most of ViTs presented in this paper
normalize image resolution by padding the patches derived from the
partitioning the image.

##  4 Methodology and Experiments

As our overarching goal is to generate the caption from the cartoon alone, we
employ multi-modal vision encoder-decoder models for this task.

As one of the main findings in [5], computer vision serves as a bottleneck for
top-quality explanation generation; models that have access to human-labeled
cartoon scene descriptions perform significantly better. However, providing
the model with both cartoon and human-annotated metadata does not accurately
mimic the task of the New Yorker Cartoon Caption Contest — the metadata for a
given cartoon (e.g. image description) is inherently up for interpretation. It
is the creative task of determining which aspects of the cartoon to
satirically emphasize that is highlighted in such a contest. Thus, providing
the model with metadata, which likely improves performance, means that we are
leaking information to the model that a normal human participant in the
contest does not have. In the remainder of this section, we outline 8
different training pipelines we devised for this task. We used the Hugging
Face and Pytorch libraries for our tasks [27, 16]

###  4.1 CLIP + GPT2

For our baseline, we implemented a CLIP-GPT2 encoder-decoder model that takes
a cartoon as input and is required to output the caption. We employ the CLIP
vision transformer model [19] to map all input cartoons into an multi-modal
embedding space. From this embedding space, we use GPT2 [20] as a decoder,
translating this image embedding into a cartoon caption by using cross-
attention between the encoder’s output and the decoder’s caption. In the
experiments described in this section, we finetune our model end-to-end on the
dataset for 10 epochs, using the AdamW optimizer with a learning rate of 5e-5.

####  4.1.1 Caption Only

In this formulation, we ignore all metadata, simply tasking the model to
generate the caption from the image alone. We include this approach as a
baseline for this task.

####  4.1.2 Caption + Metadata (all)

In this formulation, we prepend the metadata to the target label, so that the
caption is generated last. This training setup encourages the model to learn
to describe the image and generate relevant metadata first. Due to the
transformer architecture in decoder models, the model would use its attention
mechanism to refer to previously self-generated metadata, culminating in a
generated caption that is well-informed of the prior metadata. Figure 2 shows
an example of synthesized target label.

Label: |  scene: a field description: 3 dogs are out in a field. One is standing in front of the other two and preparing to throw a stick uncanny: It’s unusual to see dogs holding items as they do not have oposable digits entities: Anthropomorphism, Dog, Baseball. caption: He’s his own best friend.
---|---
Figure 2: Modified target label with metadata prepended for the cartoon in
Figure 1

###  4.2 LLaVA-NeXT

LLaVA-based models are variants of vision encoder-decoder models, which
projects the input image onto an embedding space that has the same dimension
as the language embedding space. By training this light-weight projection
layer, LLaVA-NeXT is able to turn images into tokens, which are then subject
to the same attention mechanisms as the textual data.

In [13], they used a CLIP-ViT-L/14 to Vicuna LLM encoder-decoder model.
However, in our experiments, we use the LLaVA-NeXT model. In the original
paper, the authors use a CLIP-ViT-L-336px to Vicuna LLM [12], but we make use
of a variant hosted on Hugging Face, which replaces Vicuna with Mistral LLM.

####  4.2.1 0-shot Setting

We evaluate the efficacy of LLaVA-NeXT on this task in a 0-shot setting.
Figure 3 shows the 0-shot prompt we used to query the model, and we extract
its response as the generated caption. The model is not presented any examples
in advance, and must execute the task using only its prior knowledge base.

Input: |  [INST] <image> Please write a winning caption for the provided cartoon in the New Yorker caption contest. [/INST]
---|---
Figure 3: 0-shot prompt for LLaVA-NeXT and GPT-4V. The <image> token refers to
the placement of the image embedding for the LLaVA-NeXT Model.

####  4.2.2 5-shot Setting

In addition to the 0-shot setting, we evaluate the efficacy of LLaVA-NeXT on
this task in a 5-shot setting. Specifically, we show the model 5 examples, as
demonstrated in Figure 4, with the intention to encourage the model to learn
the style of the task.

Input: |  [INST] <image1> [/INST] CAPTION 1
---|---
Input: |  [INST] <image2> [/INST] CAPTION 2
Input: |  [INST] <image3> [/INST] CAPTION 3
Input: |  [INST] <image4> [/INST] CAPTION 4
Input: |  [INST] <image5> [/INST] CAPTION 5
Input: |  [INST] <image> Please write a winning caption for the provided cartoon in the New Yorker caption contest.[/INST]
Figure 4: 5-shot prompt for LLaVA-NeXT. The <image> token refers to the
placement of the image embeddings. The final user input asks the model to
generate a caption after it has seen 5 previous examples of human-written
winning captions.

####  4.2.3 Chain-of-Thought Prompting

Chain-of-Thought (CoT) is a prompt engineering technique designed to enhance
language models’ performance on logic, calculation, and decision-making tasks
by structuring the input prompt to simulate human reasoning [26]. We provide
LLaVA-NeXT with a series of intermediate reasoning steps that encourages more
sophisticated decision-making. Figure 5 showcases our proposed Chain-of-
Thought prompting technique.

Input: |  [INST] <image> Please write an engaging, witty entry for the New Yorker Caption Contest based on this image. You are competing against 5,000 entrants from across the world, many of whom have made hundreds of previous submissions.
---|---
|  Please follow these rules:
|  1\. First, describe the comic in as much objective detail as you can,
subject to a limit of 4 sentences max. Leave no stone unturned and describe
all main objects and characters.
|  2\. Then, list as many uncanny aspects of the comic as you can think of
(list at least 5). Then, for each one, explain why it is uncanny and why the
comic artist likely inserted it. Think about what makes each aspect of the
unusual occurrence humorous.
|  3\. Based on your descriptions above and using great taste, come up with 5
potential caption concepts, each of which should be your best work and
maximally diverse. Then, roast each one in 2-3 sentences.
|  4\. Write 5 improved concepts after that, explaining why each is better.
|  [/INST]
LLM: |  GENERATED RESPONSE
Input: |  [INST] Thanks! Now, select the best caption from all the versions you’ve written above. Write only the caption as a single line, with no other text. [/INST]
Figure 5: Chain-of-Thought prompting for both LLaVA-NeXT and GPT-4V.

####  4.2.4 Finetuning with QLoRA

We finetuned LLaVA-NeXT on our training set in a question-answer format to
remain consistent with the model’s prior pre-training. Low-Rank Adaptation
(LoRA) is a popular finetuning mechanism for LLMs that significantly reduces
the number of trainable parameters without losing the benefit of downstream
finetuning [6]. LoRA approximates updates to weight matrices with a low-rank
matrix decomposition. Thus, for a given weight matrix
$W_{0}\in\mathbb{R}^{d\times k}$, LoRA constrains the update $\Delta W$ to be
approximated by a low-rank matrix decomposition $\Delta W\approx BA$, where
$B\in\mathbb{R}^{d\times r}$ and $A\in\mathbb{R}^{r\times k}$, and
$r\ll\min(d,k)$. This means that the update to weight matrices becomes

$h=(W_{0}+\Delta W)x\approx(W_{0}+BA)x.$

During training, $W_{0}$ is frozen and does not receive gradient updates,
while A and B are trainable parameters.

Due to compute limitations, we use the 7B parameter model as opposed to the
34B parameter model. In addition, we implemented quantized LoRA (QLoRA) and
LoftQ initialization [3, 8] to further optimize memory usage. We finetuned the
LLaVA-NeXT model end-to-end on an A100 80GB GPU for 10 epochs with an $r$
value of 32, learning rate of 1e-4, $\alpha$ of 16 on all linear layers.

###  4.3 GPT-4V

We repeated the same evaluation of this task using the GPT-4V multi-modal
model with both 5-shot and complex CoT mechanisms. We used the OpenAI API to
query GPT-4V and utilized the same prompts as our evaluation for LLaVA-NeXT.

GPT-4V (5-shot)
I’m here for the job interview. I have extensive experience running through
mazes. LLaVA-NeXT (0-shot)
A moment of scientific discovery, where the mouse’s curiosity meets the
scientist’s expertise. GPT-4V (complex CoT)
Should I be concerned that my new pharmacist is more interested in running
trials than running prescriptions? LLaVA-NeXT (5-shot)
As the lab assistant, I always find it strange that the mouse holds the book…
but never reads it. LLaVA-NeXT Finetuned (0-shot)
A symposium on advanced hearing aids, featuring the latest in supersonic
frequency bending technology. CLIP-GPT2 (caption only)
This is the first time I’ve had a neck injury. CLIP-GPT2 (all metadata)
Can’t you just let him finish? LLaVA-NeXT (complex CoT)
Who knew mice had such a high IQ? Figure 6: Generated Captions for New Yorker
Caption Contest #54.
Winning caption: And, when you get hungry, the cafeteria is to your right,
left, left, right, left, straight, right, straight, left, and then you push on
the big lever.

##  5 Results

We compare our generated captions to the human-written captions via automated
evaluation metrics (e.g., BLEU and ROUGE) and manual quality examination.

###  5.1 BLEU and ROUGE Scores

A well-known metric for automated quality evaluation of machine generated text
is the Bilingual Evaluation Understudy (BLEU) score metric, which is designed
to evaluate the quality of text from one natural language to another[15, 9].
All models performed poorly with respect to the BLEU score metric. Recall-
Oriented Understudy for Gisting Evaluation (ROUGE) is another popular scoring
metric for evaluating automatic summarizing and machine translation. Figure 7
shows the ROUGE scores for each of our models on the test set. While our
models performed poorly with respect to these automated metrics, we will later
discuss why these automated metrics, which are designed for machine
translation / summarization tasks, do not serve as a good metric for this
usecase.

###  5.2 Manual Quality Evaluation

To surface a more accurate metric for this task, we propose a criteria called
SS-SCORE metric. For each of the 130 test images, we compare all the machine-
generated captions with each other, using the human-written caption as a
reference, and select the best machine-generated caption. The quality of the
generated caption were based on the following criteria:

  1. 1.

Captions should contain descriptive information about the cartoon that extends
beyond listing features in the image. An observer who examines the image then
subsequently internalizes the caption should most often understand the meaning
behind why a situation was described the way it was or why a character would
have contributed their specific dialogue.

  2. 2.

Good captions do not necessarily have to exhaustively describe the situation,
but rather convey humor, wit, or a voice that aptly fits within the context of
the scene.

  3. 3.

A caption that addresses similar themes as the winning caption in the
associated contest may lead to a more favorable selection criteria on the
grounds of resembling an officially winning entry.

After tallying each test example, for each model, we aggregate the statistics
by calculating the number of times that model generated a caption that was
deemed the best, normalizing by the total number of test examples to obtain a
proportion which we equate to the SS-SCORE. These metrics can be found in
Figure 7.

##  6 Discussion

###  6.1 Automated Metrics

From Figure 7, we note that according to the various types of ROUGE scores,
the CLIP-GPT2 model performs best, and by prepending the metadata to the
target labels, we achieve a marginally higher ROUGE score. However, automated
evaluation metrics such as BLEU and ROUGE do not accurately capture the
quality of machine-generated captions. Since BLEU and ROUGE both operate on
$n$-grams, these automated metrics are capturing the commonalities between
machine generated and human-written captions at the individual word level
without regards to semantic meaning. For a task as sophisticated as this one,
word-level similarity is highly uncorrelated with caption quality because
there are multiple ways to write humorous captions. Rather than regurgitating
the content of the cartoon, contest participants are required to execute
complex reasoning that relies on social and cultural nuances related to the
cartoon, and based on the diverse range of high quality submissions for a
single image, there is no single method for generating a good caption.

###  6.2 Manual Quality Evaluation

Based on SS-SCORE, we notice that GPT-4V in the 5-shot setting performs the
best, with GPT-4V in the Chain-of-Thought setting performing second best. This
suggests that larger parameter models are able to perform better on this task,
likely because they are trained on more data and have a larger knowledge base.
GPT-4 models are constantly finetuned on the most recent data, which makes
them particularly knowledgeable about both historic and current events; this
knowledge is particularly valuable for this cartoon captioning contest, since
cartoons from The New Yorker are often timely caricatures of socio-cultural
situations.

It is particularly interesting that the 5-shot setting performs almost twice
as well as the Chain-of-Thought setting. This suggests that GPT-4V is capable
of few-shot learning, and presenting examples of cartoon-caption pairs to
GPT-4V appears more helpful for the model to adapt to the specific tone
typically found in high quality cartoon captions compared to Chain-of-Thought
prompting. Additionally, GPT-4V demonstrated other unique ways of expressing
creative generated captions such as the use of puns and double entendres (e.g.
a “full-blown service” for a tuba-playing waiter) or references to movies or
literature (e.g. “Cloudy with a Chance of Armageddon” for a cloud with a man’s
face referring to the movie Cloudy with a Chance of Meatballs).

###  6.3 Qualitative Evaluation

Figure 6 shows one particular test example with the generated captions from
all of our models. We also showcase the winning caption selected by The New
Yorker magazine.

We note that although the CLIP-GPT2 models are able to mimic the correct tone
and style of high quality captions, they are unable to capture the semantic
content of the image. This is likely because the CLIP vision model was trained
on natural images, and these cartoons do not exhibit the same color
distribution as those natural images. We found that even with fine-tuning,
CLIP is unable to capture the proper image features for GPT2 to decode.
Moreover, the CLIP-GPT2 encoder-decoder model is far from the state-of-the-art
vision transformer models, and due to the relatively small size of the model,
it is highly unlikely that it will have the vast social knowledge to create
humorous captions.

We also note that LLaVA-NeXT in the 0-shot setting has a tendency to create a
title for the cartoon rather than a caption for it. Even when we finetune
LLaVA-NeXT on the dataset, the model still has the tendency to generate a
title. We hypothesize that this is because LLaVA-NeXT does not understand the
prompt, and it has no prior knowledge of The New Yorker cartoon caption
contest. Thus, it is likely that “New Yorker” was tokenized into an unknown
token which has zero knowledge about the magazine and its history, or at the
very most has only surface-level knowledge of that linguistic entity.

However, in the 5-shot setting LLaVA-NeXT is able to generate a caption that
more closely mimics what could be found in a high quality caption. It begins
to use pronouns like in the sentence “I always find it strange…” which
indicates that it is learning to understand that high quality captions are
typically statements made by characters in the cartoon. Moreover, the Chain-
of-Thought setting shows further improvements, with the generated caption
demonstrating a clear understanding that character dialogue is important for
humorous captions.

Finally, the GPT-4V models perform the best. As discussed before, GPT-4V
models have the most recent and largest knowledge database, and are able to
successfully incorporate that knowledge into the caption. This is exhibited in
its generated captions in Figure 6, which not only shows its understanding of
the style of the caption, but also an understanding of what mice typically do
in laboratory settings. Thus, we believe that GPT-4V is the current state-of-
the-art model for this cartoon captioning task.

| ROUGE-1 | ROUGE-2 | ROUGE-L | SS-SCORE
---|---|---|---|---
GPT-4V (CoT) | 0.0653 | 0.0049 | 0.0591 | 0.19
GPT-4V (5-shot) | 0.0648 | 0.0019 | 0.0598 | 0.37
LLaVA-NeXT (0-shot) | 0.0609 | 0.0006 | 0.0568 | 0.09
LLaVA-NeXT Finetuned (0-shot) | 0.0637 | 0.0036 | 0.0575 | 0.05
LLaVA-NeXT (5-shot) | 0.0586 | 0.0024 | 0.0540 | 0.01
LLaVA-NeXT (CoT) | 0.0672 | 0.0022 | 0.0617 | 0.04
CLIP-GPT2 (caption only) | 0.0733 | 0.0062 | 0.0665 | 0.11
CLIP-GPT2 (all metadata) | 0.0735 | 0.0039 | 0.0673 | 0.14
Figure 7: Evaluation metrics for all models

##  7 Future Steps

The task of generating humorous and contextually relevant captions for New
Yorker cartoons presents unique challenges that require sophisticated AI
models. While our current work demonstrates significant advancements, several
avenues for future exploration could enhance the capabilities and performance
of such models. Our experiments with models like GPT-4V and LLaVA-NeXT
highlight the importance of model size and the breadth of pre-training data.
Future work could involve training even larger models, which inherently
possess a more extensive knowledge base and greater capacity for understanding
nuanced contexts. Utilizing models like the forthcoming GPT-5 could
potentially offer improvements in generating high-quality, humorous captions.
Additionally, leveraging distributed training on multiple GPUs or TPUs can
facilitate the training of these larger models more efficiently.

Prompt engineering plays a crucial role in guiding the model to produce
desirable outputs. Our use of few-shot learning and Chain-of-Thought (CoT)
prompting has shown promising results, but further refinement and innovation
in prompt design could yield even better performance. Especially for the
LLaVA-NeXT model, we noticed many responses misinterpreted the task at hand by
describing the image like a traditional caption more often than generating a
line that has character or wit. Another option to make the captions more
humorous and relevant is to explore hierarchical prompting, where prompts are
dynamically adjusted based on intermediate model responses.

While the New Yorker Caption Contest dataset provides a robust foundation,
curating additional datasets specifically aimed at caption generation could
offer substantial benefits. These datasets could include a broader range of
cartoons and comics, along with detailed annotations on humor, context, and
cultural references. Collaborating with cartoonists and humorists to annotate
datasets with insights into why certain captions are humorous could provide
invaluable training data for models.

## References

  * [1] A. M. Ahmed, M. Yousef, K. F. Hussain, and Y. B. Mahdy.  Enhancing image captioning with depth information using a transformer-based framework, 2023.
  * [2] A. Chandrasekaran, A. K. Vijayakumar, S. Antol, M. Bansal, D. Batra, C. L. Zitnick, and D. Parikh.  We are humor beings: Understanding and predicting visual humor, 2016.
  * [3] T. Dettmers, A. Pagnoni, A. Holtzman, and L. Zettlemoyer.  QLoRA: Efficient finetuning of quantized LLMs.  In Thirty-seventh Conference on Neural Information Processing Systems, 2023.
  * [4] A. Dosovitskiy, L. Beyer, A. Kolesnikov, D. Weissenborn, X. Zhai, T. Unterthiner, M. Dehghani, M. Minderer, G. Heigold, S. Gelly, J. Uszkoreit, and N. Houlsby.  An image is worth 16x16 words: Transformers for image recognition at scale, 2021.
  * [5] J. Hessel, A. Marasović, J. D. Hwang, L. Lee, J. Da, R. Zellers, R. Mankoff, and Y. Choi.  Do androids laugh at electric sheep? humor ”understanding” benchmarks from the new yorker caption contest, 2023.
  * [6] E. J. Hu, yelong shen, P. Wallis, Z. Allen-Zhu, Y. Li, S. Wang, L. Wang, and W. Chen.  LoRA: Low-rank adaptation of large language models.  In International Conference on Learning Representations, 2022.
  * [7] L. Jain, K. Jamieson, R. Mankoff, R. Nowak, and S. Sievert.  The New Yorker cartoon caption contest dataset, 2020.
  * [8] Y. Li, Y. Yu, C. Liang, N. Karampatziakis, P. He, W. Chen, and T. Zhao.  Loftq: LoRA-fine-tuning-aware quantization for large language models.  In The Twelfth International Conference on Learning Representations, 2024.
  * [9] C.-Y. Lin and F. J. Och.  ORANGE: a method for evaluating automatic evaluation metrics for machine translation.  In COLING 2004: Proceedings of the 20th International Conference on Computational Linguistics, pages 501–507, Geneva, Switzerland, aug 23–aug 27 2004. COLING.
  * [10] T.-Y. Lin, M. Maire, S. Belongie, L. Bourdev, R. Girshick, J. Hays, P. Perona, D. Ramanan, C. L. Zitnick, and P. Dollár.  Microsoft coco: Common objects in context, 2015.
  * [11] H. Liu, C. Li, Y. Li, and Y. J. Lee.  Improved baselines with visual instruction tuning, 2023.
  * [12] H. Liu, C. Li, Y. Li, B. Li, Y. Zhang, S. Shen, and Y. J. Lee.  Llava-next: Improved reasoning, ocr, and world knowledge, January 2024.
  * [13] H. Liu, C. Li, Q. Wu, and Y. J. Lee.  Visual instruction tuning.  ArXiv, abs/2304.08485, 2023.
  * [14] K. Nakashima, H. Kataoka, A. Matsumoto, K. Iwata, and N. Inoue.  Can vision transformers learn without natural images?, 2021.
  * [15] K. Papineni, S. Roukos, T. Ward, and W.-J. Zhu.  Bleu: a method for automatic evaluation of machine translation.  In Proceedings of the 40th Annual Meeting on Association for Computational Linguistics, ACL ’02, page 311–318, USA, 2002. Association for Computational Linguistics.
  * [16] A. Paszke, S. Gross, F. Massa, A. Lerer, J. Bradbury, G. Chanan, T. Killeen, Z. Lin, N. Gimelshein, L. Antiga, A. Desmaison, A. Köpf, E. Yang, Z. DeVito, M. Raison, A. Tejani, S. Chilamkurthy, B. Steiner, L. Fang, J. Bai, and S. Chintala.  Pytorch: An imperative style, high-performance deep learning library.  ArXiv, abs/1912.01703, 2019.
  * [17] W. Qin, R. Xu, P. Huang, X. Wu, H. Zhang, and L. Luo.  What a whole slide image can tell? subtype-guided masked transformer for pathological image captioning, 2023.
  * [18] D. Radev, A. Stent, J. Tetreault, A. Pappu, A. Iliakopoulou, A. Chanfreau, P. de Juan, J. Vallmitjana, A. Jaimes, R. Jha, and R. Mankoff.  Humor in collective discourse: Unsupervised funniness detection in the New Yorker cartoon caption contest.  In LREC, 2016.
  * [19] A. Radford, J. W. Kim, C. Hallacy, A. Ramesh, G. Goh, S. Agarwal, G. Sastry, A. Askell, P. Mishkin, J. Clark, G. Krueger, and I. Sutskever.  Learning transferable visual models from natural language supervision, 2021.
  * [20] A. Radford, J. Wu, R. Child, D. Luan, D. Amodei, and I. Sutskever.  Language models are unsupervised multitask learners.  2019\.
  * [21] D. Shahaf, E. Horvitz, and R. Mankoff.  Inside jokes: Identifying humorous cartoon captions.  In KDD, 2015.
  * [22] D. Sun, Y. Bao, and X. Cao.  A lightweight transformer for remote sensing image change captioning, 2024.
  * [23] A. Teramoto, A. Michiba, Y. Kiriyama, T. Tsukamoto, K. Imaizumi, and H. Fujita.  Automated report generation for lung cytological images using a cnn vision classifier and multiple-transformer text decoders: Preliminary study, 2024.
  * [24] A. Vaswani, N. Shazeer, N. Parmar, J. Uszkoreit, L. Jones, A. N. Gomez, L. Kaiser, and I. Polosukhin.  Attention is all you need, 2023.
  * [25] R. Vedantam, X. Lin, T. Batra, C. L. Zitnick, and D. Parikh.  Learning common sense through visual abstraction.  In 2015 IEEE International Conference on Computer Vision (ICCV), pages 2542–2550, 2015.
  * [26] J. Wei, X. Wang, D. Schuurmans, M. Bosma, b. ichter, F. Xia, E. Chi, Q. V. Le, and D. Zhou.  Chain-of-thought prompting elicits reasoning in large language models.  In S. Koyejo, S. Mohamed, A. Agarwal, D. Belgrave, K. Cho, and A. Oh, editors, Advances in Neural Information Processing Systems, volume 35, pages 24824–24837. Curran Associates, Inc., 2022.
  * [27] T. Wolf, L. Debut, V. Sanh, J. Chaumond, C. Delangue, A. Moi, P. Cistac, T. Rault, R. Louf, M. Funtowicz, J. Davison, S. Shleifer, P. von Platen, C. Ma, Y. Jernite, J. Plu, C. Xu, T. Le Scao, S. Gugger, M. Drame, Q. Lhoest, and A. Rush.  Transformers: State-of-the-art natural language processing.  In Q. Liu and D. Schlangen, editors, Proceedings of the 2020 Conference on Empirical Methods in Natural Language Processing: System Demonstrations, pages 38–45, Online, Oct. 2020. Association for Computational Linguistics.
  * [28] K. Xu, J. Ba, R. Kiros, K. Cho, A. Courville, R. Salakhutdinov, R. Zemel, and Y. Bengio.  Show, attend and tell: Neural image caption generation with visual attention, 2016.
  * [29] C. L. Zitnick, D. Parikh, and L. Vanderwende.  Learning the visual interpretation of sentences.  In 2013 IEEE International Conference on Computer Vision, pages 1681–1688, 2013.