# Unsupervised authorship attribution

David Fifield, Torbjørn Follan, Emil Lunde

We describe a technique for attributing parts of a written text to a set of
unknown authors. Nothing is assumed to be known a priori about the writing
styles of potential authors. We use multiple independent clusterings of an
input text to identify parts that are similar and dissimilar to one another.
We describe algorithms necessary to combine the multiple clusterings into a
meaningful output. We show results of the application of the technique on
texts having multiple writing styles.

##  1 Problem statement

Accept as input a written text and a number $n$ of presumed authors. For each
small segment of the text—for example, for each word or line—output a weight
assignment for the labels $1,\ldots,n$, with the weights summing to $1$.

Nothing is previously assumed about the authors—there are no training samples
of writing by potential authors. The author labels in the output are just
integers, not tied to any real personal identity. The weights on each segment
may be interpreted as a probability distribution over authors; for example a
segment may appear as 15% Unknown Author 1, 80% Unknown Author 2, and 5%
Unknown Author 3. The output labels are arbitrary: any permutation of labels
that is applied consistently across all segments is as good as any other.

The input text is an ordinary written document in a human language such as
English. The text doesn’t have explicit section markings or other divisions;
from the point of view of the algorithm it is just a long sequence of words.
The $n$ authors may have worked together on the document at the same time, or
the text may have been assembled using parts from different authors separated
in time.

###  1.1 Prior work

There is a great deal of literature on the topic of supervised authorship
attribution—extracting lexical and other features from writing samples known
to belong to different authors, and using those to help decide the authorship
of an unknown text. Stamatatos [13] gives a survey. There is comparatively
little on assigning authorship to authors without prior knowledge of their
writing styles. Koppel et al. [11] use synonym choice, among other features,
to learn authorship divisions in input texts artificially constructed from
different parts of the Bible.

##  2 Solution technique

A technique for assigning $n$ labels to parts of a document suggests itself:
Break the input into fragments, cluster the fragments into $n$ clusters, and
assign a label to each fragment equal to the cluster label. The clustering may
be based on, for example, stylometric features extracted from the fragments.
We expect that fragments within a cluster are similar to one another but
dissimilar to fragments in other clusters. If the similarities and
dissimilarities come principally from different authors’ stylistic choices,
the cluster labels will reflect authorship.

There is a problem with the technique just described: A fragment of text needs
to have a certain minimum length in order to represent an author’s style. The
minimum has been found to vary with different classification methods [13], but
it should be at least several hundred words for English text. Longer fragments
convey more of an author’s style, at the cost of low-resolution output: If all
fragments are 1000 words long, the output may show a change in authorship only
every 1000 words. Smaller fragments give better resolution, at the cost of
less distinctive style within fragments and less robust clustering.

Our proposed solution to this problem is to repeat the naive clustering $m$
times, each time using fragments that are shifted relative to the previous
time. Reassign cluster labels within clusterings to make the clusterings agree
as much as possible. Output the average of the $m$ relabeled clusterings. This
is a summary of the procedure:

  1. 1.

Make $m$ copies of the input, and break each copy into fragments, with the
fragments in one copy being shifted relative to the other copies.

  2. 2.

Cluster each set of fragments independently according to some similarity
measure.

  3. 3.

Permute the labels within each clustering to maximize pairwise clustering
agreement.

  4. 4.

Output the average of the relabeled clusterings.

For example, let the fragment size be 1000 words and the shift offset be 50
words. These parameters lead to $m=1000/50=20$ sets of fragments. The first
comprises words 1–1000,1001–2000,…; the second words 51–1050,1051–2050,…; and
so on. We omit short fragments, for example the one composed of words 1–51
that would otherwise be in the second fragment set. The output will be short
segments of only 50 words, each with its own weight distribution.

Fragments do not necessarily have to be the same size. They should be long
enough to represent a single author’s style, but not so long that fragments
tend to have more than one or two authors. Fragments that straddle a true
change in authorship tend to be assigned to one or the other author according
to how much each author contributed to the fragment. The act of averaging
fragments that start at different offsets therefore tends to cause label
weights to transition smoothly between one label and another, at least where
the change in authorship is clear.

The overall effect is that we get output fine-grained enough to localize
authorship changes, while still doing clustering over longer, stylistically
meaningful fragments. Repeated clustering introduces an additional
complication, taken care of in step 3. Because the $m$ clusterings are
independent, they will not in general use the same label for the same unknown
author. Step 3 permutes labels within each clustering, as we are free to do,
to make the labels agree as much as possible.

The following sections consider different techniques for doing step 2
(clustering), and step 3 (matching of clustering labels).

##  3 Clustering techniques

In this section we consider ways to cluster text fragments according to their
stylistic similarity. Clustering over feature vectors (Section 3.1) is
amenable to a variety of machine learning techniques such as dimensionality
reduction (Section 3.2); it requires some amount of natural language
processing in order to extract the features. Compression clustering (Section
3.3) does not deal with meaningful feature vectors, but may be used even when
it is inconvenient to extract features, for example when the input is in a
language for which natural language processing libraries are lacking.

###  3.1 Stylometric feature vectors

There are many ways to transform text into a numeric feature vector. Typical
features are things like average word length, and frequency of use of
punctuation—but there exist much more sophisticated features. Once a set of
feature vectors has been produced, they may be clustered using any of a
variety of general-purpose clustering algorithms.

JStylo [12] is a framework for authorship attribution, usually used to
evaluate a test set of documents against a training set. We used JStylo to
generate feature vectors over fragment sets that were previously generated. We
used the built-in “WritePrints (Limited)” feature set. The WritePrints feature
set is language-aware and specific to English, having features not only for
things like word length and character frequency, but also parts of speech and
frequency of function words [2]. We then used the machine learning framework
Weka [9] to generate $m$ sets of $n$ clusters using the EM clusterer.
Clustering of JStylo-produced feature vectors was used to get the results in
Section 5.1.

###  3.2 Dimensionality reduction

Dimensionality reduction in the sense of clustering, is the act of reducing
the potentially large number of features in a feature vector. This is used to
speed up the clustering software. The algorithm we use is a random-projection-
type $k$-means algorithm proven to distort the clusters by at most a factor of
$2+\epsilon$ for an error parameter $\epsilon$. The algorithm we use is due to
Boutsidis, Drineas, and Zouzias [6, Algorithm 1]. It uses a random scaled $\pm
1$ $d\times t$ matrix to reduce $d$ dimensions to $t$ dimensions, then runs an
approximate $k$-means algorithm on the projected vectors.

The parameter $t$, which is the new number of features, is independent of $d$,
which is the original number of features. In our case, $d$, the number of
features, was around 450–500, and $k$, the number of clusters, was $2$. We
chose $\epsilon=0.2$ and $c=2$, giving $t=100$ output features.

###  3.3 Compression clustering

Clustering based on traditional stylometric features requires knowledge about
the data in order to extract useful feature that can be used to compare the
objects. This is a problem since there are a lot of different languages, with
different properties and alphabets, that we may want to cluster.

Clustering by compression is a clustering method introduced by Cilibrasi and
Vitanyi [7] that doesn’t require any subject-specific features or knowledge
about the data. The idea is that two objects are deemed similar if one can use
information about one object to compress the other. More formally the
normalized compression distance between objects $x$ and $y$ is

$\mathrm{NCD}(x,y)=\frac{\mathrm{C}(xy)-\min(\mathrm{C}(x),\mathrm{C}(y))}{\max(\mathrm{C}(x),\mathrm{C}(y)))},$

where $\mathrm{C}(x)$ is the compressed size of object $x$, and $xy$ is the
concatenation of $x$ and $y$.

The normalized compression distance allows us to cluster texts written in
languages we don’t have tools to extract features from. To test this method we
used the CompLearn utilities [1] to generate a distance matrix between
fragments and the $k$-medioids algorithm to cluster them. An example of the
application of this approach appears in Section 5.2.

##  4 Clustering matching algorithms

The output of step 2 of the algorithm is a set of $n$ clusterings, each one
being a labeling of fragments. The labels within each clustering may be
permuted arbitrarily. The job of step 3 is to find a permutation of labels
within each clustering so that the labels across clusterings agree as much as
possible.

We examined two special cases of this optimization: the case where there are
two authors and multiple clusterings, and the case where there are multiple
authors and two clusterings.

###  4.1 Agreement between labeled fragments

We define the normalized _agreement_ between a pair of fragment sets $A$ and
$B$ on a text $T$:

$\mathrm{agreement}(A,B)=\frac{1}{|T|}\sum_{f_{a}\in A}\sum_{f_{b}\in
B}|f_{a}\cap f_{b}|\sum_{i=1}^{n}(f_{a})_{i}(f_{b})_{i}.$

$|T|$ is the length of the text, for example a count of its words. $A$ and $B$
are fragment sets; $f_{a}$ and $f_{b}$ are individual fragments within them.
$|f_{a}\cap f_{b}|$ is the size of the intersection of $f_{a}$ and $f_{b}$;
for nonoverlapping fragments it is $0$. $(f_{A})_{i}$ and $(f_{B})_{i}$ are
the weight assignments for author $i$ in $f_{A}$ and $f_{B}$ respectively.
Simply stated, the agreement between two fragments is the size of their
intersection times the dot product of their weight vectors. The total
agreement between two sets is the sum of agreements between all fragments in
$A$ and all fragments in $B$; it is a number between $0$ and $1$.

The sum-of-pairs agreement of a sequence $M$ containing $m$ sets of fragments
is the normalized sum of their pairwise agreements:

$\mathrm{agreement}(M)=\frac{1}{{m\choose 2}}\sum_{1\leq i<m}\sum_{i<j\leq
m}\mathrm{agreement}(M_{i},M_{j}).$

The use of a sum-of-pairs objective function has analogs in algorithms used in
biology for multiple sequence alignment—there, the algorithm must choose where
to insert gaps into various sequences in order to maximize a score; here, the
algorithm’s choice is over permutations of labels.

We use the agreement objective for two purposes. First, sum-of-pairs agreement
is used to optimize labels of the $m$ independent fragment clusterings.
Second, single-pair agreement is used to assess the effectiveness of the
algorithm’s output in comparison with some known ground truth. Examples of the
latter use are given in Section 5. In the specific case where there are two
authors, and label weights are constrained to be $0$ or $1$, optimization of
single-pair agreement has a natural interpretation as a binary classification
problem. Without loss of generality, fix the labels of one fragment set in the
pair. The fragments with label $1$ are the positive instances, and those with
label $2$ are the negative instances. Then the agreement is just the
normalized sum of true positives and true negatives over the total:
$\frac{1}{|T|}(\mathrm{TP}+\mathrm{TN})/(\mathrm{TP}+\mathrm{TN}+\mathrm{FP}+\mathrm{FN})$.

###  4.2 Two clusterings, $n>2$ authors: maximum-weight matching

When there are only two clusterings, the problem of how to permute the two
label sets for maximum agreement is an instance of maximum-weight matching.
Build a complete bipartite graph with $n$ vertices on each side, corresponding
to the $n$ labels in each clustering. For each pair $1\leq i,j\leq n$, assign
a weight equal to the size of the intersection between fragments labeled $i$
on the left and fragments labeled $j$ on the right; this is the size of the
contribution toward the total agreement that would be achieved if label $i$ on
the left were changed to $j$. There are algorithms to solve maximum-weight
matching optimally in polynomial time.

###  4.3 $m>2$ clusterings, two authors: MAX-CUT using a semidefinite
programming relaxation

There are only two permutations of two labels: $(1,2)$ and $(2,1)$. The choice
of a permutation of labels in a fragment set therefore boils down to the
choice of whether to flip the labels, or not. Consider a two-element weight
vector within a fragment $(w_{1},w_{2})$. Because $w_{1}+w_{2}=1$, it is the
case that flipping the labels has the same effect as flipping the weights: the
weight vector becomes $((1-w_{1}),(1-w_{2}))=(w_{2},w_{1})$. The same is true
of the agreement, which is the sum of dot products of weight vectors:

$\mathrm{agreement}(1-A,B)=\mathrm{agreement}(A,1-B)=1-\mathrm{agreement}(A,B).$

We seek to decide, for each fragment set, whether to flip its labels or not,
so that the sum or pairwise agreements is as large as possible.

Optimizing the agreement between multiple two-author clusterings can be posed
as an instance of the MAX-CUT problem. An instance of MAX-CUT is an undirected
graph $G=(V,E)$ with nonnegative edge weights. A solution to the problem is a
mapping from vertices to labels in $\\{-1,+1\\}$ such that the sum of weights
on edges whose endpoints have different labels is maximized. Vertices with
label $-1$ correspond to one side of a graph cut; those with label $+1$
correspond to the other side, and the objective is to maximize the weight of
the edges crossing the cut. Let an instance of MAX-CUT be given with vertex
set $V=\\{v_{1},\ldots,v_{m}\\}$ and edge weights $w_{ij}=w_{ji}$. Its
solution is the solution to the integer program:

$\displaystyle\max\quad$ $\displaystyle\sum_{0\leq i<n}\sum_{i<j\leq
m}w_{ij}\frac{1-y_{i}y_{j}}{2}$
$\displaystyle\mathrm{subject\leavevmode\nobreak\ to\quad}$ $\displaystyle
y_{i}\in\\{-1,+1\\}\qquad\forall i\in V$

The objective function has the same solution—though not the same value—after
removing constant positive factors and additive terms: $\max\sum_{0\leq
i<n}\sum_{i<j\leq m}-w_{ij}y_{i}y_{j}$ (using payoff multipliers in
$\\{-1,+1\\}$ rather than $\\{0,1\\}$).

The reduction from two-author clustering matching to MAX-CUT is achieved
simply by setting $w_{ij}=-\mathrm{agreement}(M_{i},M_{j})$. Intuitively, if
the agreement is large (greater than $0.5$), then we want $-w_{ij}y_{i}y_{j}$
to be small, meaning that $y_{i}$ and $y_{j}$ should have the same sign (flip
neither of $v_{i}$ and $v_{j}$, or flip both). If the agreement is small (less
than $0.5$), $y_{i}$ and $y_{j}$ should have opposite signs (flip one label
set of the other, but not both). More formally, we may work out a payoff
function that pays $+\mathrm{agreement}(M_{i},M_{j})$ when $y_{i}$ and $y_{j}$
have the same sign, and $-\mathrm{agreement}(M_{i},M_{j})$ when $y_{i}$ and
$y_{j}$ have opposite signs:

$\displaystyle\mathrm{agreement}(M_{i},M_{j})\left(\frac{1+y_{i}y_{j}}{2}\right)-\mathrm{agreement}(M_{i},M_{j})\left(1-\frac{1+y_{i}y_{j}}{2}\right)$
$\displaystyle=\leavevmode\nobreak\ $
$\displaystyle\mathrm{agreement}(M_{i},M_{j})\left(\frac{1+y_{i}y_{j}}{2}\right)-\mathrm{agreement}(M_{i},M_{j})\left(\frac{1-y_{i}y_{j}}{2}\right)$
$\displaystyle=\leavevmode\nobreak\ $
$\displaystyle\mathrm{agreement}(M_{i},M_{j})\left(\frac{1+y_{i}y_{j}-1+y_{i}y_{j}}{2}\right)$
$\displaystyle=\leavevmode\nobreak\ $
$\displaystyle\mathrm{agreement}(M_{i},M_{j})y_{i}y_{j}$

which has the same maximum as the MAX-CUT program when assigning
$w_{ij}=\mathrm{agreement}(M_{i},M_{j})$.

Goemans and Williamson [8] show a randomized polynomial-time
$0.878$-approximate algorithm to solve MAX-CUT by way of a semidefinite
programming relaxation. We apply this algorithm after reducing two-author
clustering matching to MAX-CUT. The output labels
$y_{1},\ldots,y_{m}\in\\{-1,+1\\}$ map to different sides of the cut. The
clusterings on one side of the cut have their labels flipped; and those on the
other side are left unchanged.

There are other natural interpretations of the edge weights we assign: If we
regard the labels as being in $\\{0,1\\}$, with aligned fragments of identical
length, then the weight is the negative of the Hamming distance between the
bit strings representing the labels. If we regard the labels as being in
$\\{-1,+1\\}$, under the same conditions, then the weight is the negative of
the dot product or cosine similarity between the label vectors.

###  4.4 $m>2$ clusterings, $n>2$ authors

We did not try to solve any matchings where there are more than two authors
and more than two clusterings. It is likely that approximate or heuristic
techniques will be necessary. We believe that the application of algorithms
used for sequence alignment in biology is promising.

##  5 Results

###  5.1 Ulysses

> Ineluctable modality of the visible: at least that if no more, thought
> through my eyes. Signatures of all things I am here to read, seaspawn and
> seawrack, the nearing tide, that rusty boot. Beginning of episode 3

The novel Ulysses by James Joyce [10] is known for the drastic and deliberate
stylistic shifts between the 18 chapters or “episodes” of which it is
composed. We used the text of the novel as an input, treating each episode as
if it had been composed by a different author. We wanted to see whether these
stylistic divisions could be automatically detected by a program.

Our input was Project Gutenberg’s copy of a pre-1923 edition of the novel,
released in July 2003. After editing to remove some episode markers, the text
has 32309 words. The shortest episode is number 2, with 4414 words, and the
longest is number 15, with 38161.

Figure 1: Attribution of Ulysses. The agreement between the computed and true
labels is $0.55$.

Figure 1 shows the process of attribution. We followed the “two clusterings,
$n>2$ authors” plan from Section 4.2. We broke the input into uniform
fragments of 5000 words, with one set of fragments offset by 2500 words
relative to the other. The output is therefore a distribution over author
labels for every 2500-word segment, where the weights on each segment are
either fully on one author or evenly split between two authors. The two
fragment sets are represented by the two gray columns at the left of the
figure. We then generated feature vectors for each 5000-word fragment using
JStylo, and clustered both fragment sets independently; this is the second
pair of columns in the figure. Author labels are shown as different colors.
Notice how there is hardly any agreement between the two columns. We then
computed an maximum-weight matching between the two fragment sets. As shown in
the third pair of columns, this makes many more labels agree. On the right of
the figure is the mean of the two relabeled clusterings, aside the known true
episode division.

We see that some episodes, perhaps those most stylistically distinct, are very
closely matched: 8, 9, 10, 15, 16, and 17. In other places the algorithm
appears confused; reasons why may include stylistic similarity between
episodes, short episodes, or episode boundaries that do not fall close to a
multiple of 2500 words. The agreement between the computed authorship and the
actual episode division is $0.55$ (random assignment would have $1/18\approx
0.06$).

###  5.2 The Homeric Hymn to Apollo

>
>
>
>
> µµ µ ,
>  µ µ :
>  µ
>  , µ .
>
> Λινες 1–4 οφ τηε ηψμν
>
>
>
> Τηε Ηομερις Ηψμνς αρε α ςολλεςτιον οφ ανςιεντ Γρεεϰ ηψμνς ηονορινγ τηε γοδς.
> Ονε οφ τηεςε, τηε Ηψμν το Απολλο, βελιεεδ το ηαε βεεν ωριττεν βεφορε τηε 5τη
> ςεντυρψ Β.῝.,
> ις υςυαλλψ πρεςεντεδ ας ονε ποεμ ιν τωο παρτς: α ςηορτερ ``Δελιαν´´ παρτ ανδ
> α λονγερ ``Πψτηιαν´´ παρτ.
> Τεξτυαλ ανδ ηιςτοριςαλ ειδενςε ςυγγεςτς τηατ, ωηιλε τραδιτιοναλλψ ςομβινεδ,
> τηεψ ωερε οριγιναλλψ ςομποςεδ ινδεπενδεντλψ.
> Αλλεν ανδ Σιϰες [3] ωριτε:
>
>
>

>>
>
>>

>> Τηε ηψμν το Απολλο, ιν ιτς πρεςεντ φορμ, μαψ βε ρεαδ ας α ςοντινυους ποεμ.
> Βυτ τηε ςοντινυιτψ λιες ονλψ ον τηε ςυρφαςε, ανδ εεν τηε μοςτ ςαςυαλ ρεαδερ
> ςαννοτ φαιλ το βε ςτρυςϰ βψ
> τηε αβρυπτ τρανςιτιον ατ . 179… Τηε υνιτψ οφ τηε ηψμν ηας βεεν δενιεδ ον
> αρτιςτις ανδ λιτεραρψ
> γρουνδς. … Ιτ τηερεφορε φολλοως τηατ τηε ηψμν ις α ςομπιλατιον οφ
> _ατ λεαςτ_ τωο οριγιναλλψ ινδεπενδεντ ποεμς.
>
>
>
>
>
>
>
>
>
> Ουρ γοαλ ωας το ςεε ιφ α ςομπυτερ ςαν αυτοματιςαλλψ δετεςτ ςτψλιςτις ςηανγε
> βετωεεν τηε τωο παρτς,
> ωιτηουτ α πριορ ςαμπλες οφ ωριτινγ φρομ πρεςυμεδ αυτηορς.
> Ωε υςεδ ας ινπυτ τηε ϒΤΦ-8 Γρεεϰ τεξτ οφ τηε ηψμν [4]
> φρομ τηε Περςευς Διγιταλ Λιβραρψ.
> Ινςλυδινγ βραςϰετεδ λινες, τηε Γρεεϰ ηψμν ςονςιςτς οφ 3840 ωορδς ιν 552
> λινες,
> ωιτη 1241 ωορδς ιν τηε Δελιαν πορτιον ανδ 2599 ιν τηε Πψτηιαν.
>
>
>
>  Φιγυρε 1: Αττριβυτιον οφ τηε Ηψμν το Απολλο.
> Τηε αγρεεμεντ βετωεεν τηε ςομπυτεδ ανδ τρυε λαβελς ις $0.84$.
>
>
>
> Φιγυρε 5.2 ιλλυςτρατες τηε προςεςς οφ ςομπυτινγ αυτηορςηιπ.
> Ωε βεγιν βψ βρεαϰινγ τηε ινπυτ ιντο φραγμεντς.
> Τηε ινπυτ ις ρατηερ ςηορτ ας αυτηορςηιπ αττριβυτιον γοες—ονλψ 3736 ωορδς.
> Ωε υςε φραγμεντς οφ 900 ωορδς, ανδ ομιτ ςηορτ φραγμεντς ατ τηε βεγιννινγ ανδ
> ενδ,
> μεανινγ τηατ εαςη ςλυςτερινγ ωιλλ ςλυςτερ ονλψ τηρεε ορ φουρ φραγμεντς.
> Τηεψ ςταρτ ατ οφφςετς οφ 90 ωορδς ρελατιε το εαςη οτηερ, ςο τηερε αρε 10
> ςλυςτερινγς.
> Τηε φιρςτ ςλυςτερινγ ςταρτς ατ οφφςετ $0$, οερ τηε φουρ φραγμεντς ςονςιςτινγ
> οφ ωορδς
> 1–900, 901–1800, 1801–2700, ανδ 2701–3600.
> Τηε νεξτ ςλυςτερινγ ις οερ ωορδς 91–990, 991–1890, 1891–2790, 2791–3690, ανδ
> ςο ον.
> Ιν τηε φιγυρε, εαςη ςολυμν ις αν ινδεπενδεντλψ ςομπυτεδ ςλυςτερινγ.
>
> Τηε ςλυςτερινγ μαψ ςομε ουτ ας ιν τηε λεφτ ςιδε οφ Φιγυρε 5.2.
> Εαςη φραγμεντ ις ςολορεδ αςςορδινγ το ιτς λαβελ αςςιγνμεντ.
> Βεςαυςε τηε αλγοριτημ ις υνςυπεριςεδ, α λιγητ βλυε λαβελ ςτανδς φορ ϒνϰνοων
> Αυτηορ 1
> ανδ α δαρϰ βλυε λαβελ φορ ϒνϰνοων Αυτηορ 2. Βεςαυςε τηε ςλυςτερινγς αρε
> ινδεπενδεντ,
> τηερε ις νοτ μυςη αγρεεμεντ βετωεεν τηεμ.
> Τηε μιδδλε παρτ οφ τηε φιγυρε ςηοως τηε ςλυςτερινγς αφτερ ρελαβελινγ βψ τηε
> ΜΑΞ-῝ϒΤ αππροξιματιον.
> Φορ εαςη ςλυςτερινγ, τηε ΜΑΞ-῝ϒΤ αππροξιματιον ηας φρεεδομ ονλψ το ινερτ τηε
> λαβελς, ορ νοτ.
> Ηερε τηε αλγοριτημ ηας ςηοςεν το ινερτ τηε λαβελς οφ τηε τηιρδ, ςιξτη,
> ςεεντη, νιντη, ανδ τεντη ςλυςτερινγς φορ βεττερ αγρεεμεντ.
> (Ινερτινγ τηε φιρςτ, ςεςονδ, φουρτη, φιφτη, ανδ ειγητη ινςτεαδ ωουλδ ρεςυλτ
> ιν τηε ςαμε ςυτ ανδ βε ϑυςτ ας γοοδ.)
>
> Τηε ριγητ ςιδε οφ Φιγυρε 5.2 ςηοως τηε μεαν οφ τηε ρελαβελεδ ςλυςτερινγ.
> Τηε βαλανςε οφ αυτηορςηιπ μαψ ςηανγε ας οφτεν ας εερψ 90 ωορδς,
> μαϰινγ ποςςιβλε τηε ςμοοτη γραδιεντ-λιϰε τρανςιτιον ςεεν ηερε.
> ῝ομπαρινγ τηε ςομπυτεδ αυτηορςηιπ ωιτη τηε ϰνοων διιςιον ιντο Δελιαν ανδ
> Πψτηιαν παρτς,
> ωε φινδ τηεμ μοςτλψ ιν αγρεεμεντ, ςταρτινγ ωιτη ονε αυτηορ ατ τηε βεγιννινγ
> ανδ τρανςιτιονινγ το ανοτηερ αβουτ 30% τηρουγη,
> τηουγη ιτ ρεμαινς ονλψ 80% ςυρε οφ τηε ςεςονδ αυτηορ υντιλ νεαρ τηε ενδ οφ
> τηε τεξτ.
> Τηε αγρεεμεντ βετωεεν τηε ςομπυτεδ αυτηορςηιπ ανδ τηε αςτυαλ Δελιαν/Πψτηιαν
> ςπλιτ ις $0.84$
> (ρανδομ αςςιγνμεντ ωουλδ ηαε $0.5$).
>
>
>
>  Φιγυρε 1: Αττριβυτιον οφ α μανγλεδ Ηψμν το Απολλο.
> Τηε Δελιαν πορτιον ηας βεεν αρβιτραριλψ παςτεδ ιντο τηε μιδδλε οφ τηε
> Πψτηιαν.
> Τηε αγρεεμεντ ις $0.72$.
>
>
>
> Τηε μεαν αςςιγνμεντ ςαννοτ ηαε α ςηαρπ βρεαϰ βετωεεν αυτηορς ας τηε τρυε
> αυτηορςηιπ δοες·
> ιτ αλωαψς ταϰες ατ λεαςτ 900 ωορδς το τρανςιτιον φρομ ονε το ανοτηερ, δυρινγ
> ωηιςη τρανςιτιον ιτ ις ον αεραγε ατ λεαςτ ηαλφ ωρονγ.
> ῝ορρεςπονδινγλψ αν υππερ βουνδ ον τηε αγρεεμεντ αςηιεαβλε υςινγ τηις
> τεςηνιχυε, ιν τηις ςαςε, ις $(3726-450)/3726\approx 0.88$.
>
> Ωε διδ ανοτηερ, αρτιφιςιαλ, τεςτ ον τηε ηψμν.
> Ιν ορδερ το προδυςε α τρυε αςςιγνμεντ τηατ ις νον-ςοντιγυους,
> ωε ςυτ τηε Δελιαν παρτ ανδ παςτεδ ιτ ιντο τηε μιδδλε οφ τηε Πψτηιαν παρτ.
> Τηε προςεςς ις ςηοων ιν Φιγυρε 5.2.
> Τηε ςεςονδ αυτηορςηιπ ςηανγε ις ερψ ςηαρπ βυτ τηε φιρςτ ις λεςς ςλεαρ.
> Τηε τοταλ αγρεεμεντ ωιτη τηε ϰνοων λαβελινγ ις $0.72$.
> Βεςαυςε τηερε αρε τωο πλαςες ωηερε αυτηορςηιπ ςηανγες, ςεπαρατεδ βψ
> μορε τηαν 900 ωορδς, τηε μαξιμυμ αγρεεμεντ ατταιναβλε βψ τηε αλγοριτημ ις
> $(3726-900)/3726\approx 0.76$.
>
>
>
> ##  Διμενςιοναλιτψ ρεδυςτιον
>
>
>
>
>
>
> Ιν ορδερ το εξπεριμεντ ωιτη διμενςιοναλιτψ ρεδυςτιον ον φεατυρε εςτορς,
> ωε νεεδ νυμερις φεατυρε εςτορς ςυςη ας τηοςε φρομ Σεςτιον 3.1.
> Ωε υςεδ ΘΣτψλο το προδυςε φεατυρε εςτορς φρομ αν Ενγλιςη τρανςλατιον οφ τηε
> Ηψμν το Απολλο [5].
> Τηις τεξτ ηας 1768 ωορδς ιν τηε Δελιαν πορτιον ανδ 3726 ιν τηε Πψτηιαν.
>
>
>
>
>
>
>
>
>
>
> Φιγυρε 1: Βεφορε ανδ αφτερ διμενςιοναλιτψ ρεδυςτιον οφ αν Ενγλιςη
> τρανςλατιον οφ τηε Ηψμν το Απολλο βψ ρανδομ προϑεςτιον.
> Τηε αγρεεμεντ βετωεεν τηε ουτπυτ λαβελς ις 0.56.
>
>
>
>
>
>
>
> Φιγυρε 5.2 ιλλυςτρατες τηε διφφερενςε ωιτη ανδ ωιτηουτ διμενςιοναλιτψ
> ρεδυςτιον.
> Νοτιςε τηατ τηεψ αγρεε μοςτλψ ιν τηε ςταρτ ανδ ατ τηε ενδ, βυτ τηεψ αρψ α
> βιτ μορε ιν τηε μιδδλε. Τηις ις πρεςυμαβλψ βεςαυςε τηε διμενςιοναλιτψ
> ρεδυςτιον ις αν $2+\epsilon$ αππροξιματιον.
>
>
>
> ##  Σουρςε ςοδε
>
>
>
>
>
> Ουρ ςουρςε ςοδε ανδ δατα φιλες αρε ααιλαβλε ιν α Γιτ ρεποςιτορψ ατ
> ηττπς://ρεπο.εεςς.βερϰελεψ.εδυ/γιτ-ανον/υςερς/φιφιελδ/ςς270.γιτ.
>
>
>
> ##  Φυτυρε ωορϰ
>
>
>
>
>
> Ωε ςηοωεδ τηατ τηε προβλεμ οφ ματςηινγ μυλτιπλε τωο-αυτηορ ςλυςτερινγς
> ρεδυςες το ΜΑΞ-῝ϒΤ,
> βυτ διδν´τ ςηοω τηατ τηε προβλεμ ις νεςεςςαριλψ ας ηαρδ ας ΜΑΞ-῝ϒΤ.
> Τηε εδγε ωειγητς ιν τηε ρεδυςτιον ηαε ςονςτραιντς· φορ εξαμπλε, τηεψ οβεψ
> τηε τριανγλε ινεχυαλιτψ.
> Ωε ωουλδ λιϰε το δετερμινε ωηετηερ τηε προβλεμ ςαν βε ςολεδ ας ωελλ ορ
> βεττερ ωιτηουτ ρεδυςινγ το ΜΑΞ-῝ϒΤ.
>
> Τηε ςαςε οφ $m>2$ ςλυςτερινγς ανδ $n>2$ αυτηορς νεεδς ινεςτιγατιον το ςεε
> ηοω ωελλ ιτ ςαν βε αππροξιματεδ.
>
>
>
> ##  Αςϰνοωλεδγεμεντ
>
>
>
>
>
> Ωε τηανϰ Θαμες ῝οοϰ ανδ Στεπηεν Σανςομ φορ τηειρ ςομμεντς ον δραφτς οφ τηις
> ρεπορτ.
>
>
>
> ## Ρεφερενςες
>
>   *
>
>
>   * [1]
> ῝ομπΛεαρν ηομε παγε.
>  ηττπ://ωωω.ςομπλεαρν.οργ/.
>
>
>   * [2]
> Αημεδ Αββαςι ανδ Ηςινςηυν ῝ηεν.
>  Ωριτεπριντς: Α ςτψλομετρις αππροαςη το ιδεντιτψ-λεελ ιδεντιφιςατιον
> ανδ ςιμιλαριτψ δετεςτιον ιν ςψβερςπαςε.
>  Α῝Μ Τρανς. Ινφ. Σψςτ., 26(2):7:1–7:29, Απριλ 2008.
>
>
>   * [3]
> Τηομας Ω. Αλλεν ανδ Ε. Ε. Σιϰες.
>  Τηε Ηομερις Ηψμνς.
>  1904\.
>
>
> ηττπ://ωωω.περςευς.τυφτς.εδυ/ηοππερ/τεξτ;δος=Περςευς:τεξτ:1999.04.0029:τεξτ=ςομμ:ποεμ=3.
>
>
>   * [4]
> Ανονψμους.
>  Ηψμν το Απολλο.
>
>
> ηττπ://ωωω.περςευς.τυφτς.εδυ/ηοππερ/τεξτ;δος=Περςευς:τεξτ:1999.01.0137:ηψμν=3.
>  Αςςεςςεδ Νοεμβερ 13, 2013.
>
>
>   * [5]
> Ανονψμους.
>  Ηψμν το Απολλο.
>
>
> ηττπ://ωωω.περςευς.τυφτς.εδυ/ηοππερ/τεξτ;δος=Περςευς:τεξτ:1999.01.0138:ηψμν=3.
>  Τρανςλατεδ βψ Ηυγη Γ. Εελψν-Ωηιτε. Αςςεςςεδ Νοεμβερ 14, 2013.
>
>
>   * [6]
> ῝ηριςτος Βουτςιδις, Αναςταςιος Ζουζιας, ανδ Πετρος Δρινεας.
>  Ρανδομ προϑεςτιονς φορ $k$-μεανς ςλυςτερινγ.
>  ῝οΡΡ, αβς/1011.4632, 2010.
>  ηττπ://βοοϰς.νιπς.ςς/παπερς/φιλες/νιπς23/ΝΙΠΣ2010῍0113.πδφ.
>
>
>   * [7]
> Ρ. ῝ιλιβραςι ανδ Π.Μ.Β. ῞ιτανψι.
>  ῝λυςτερινγ βψ ςομπρεςςιον.
>  Ινφορματιον Τηεορψ, ΙΕΕΕ Τρανςαςτιονς ον, 51(4):1523–1545,
> 2005\.
>  ηττπ://ηομεπαγες.ςωι.νλ/῀παυλ/παπερς/ςλυςτερ.πδφ.
>
>
>   * [8]
> Μιςηελ Ξ. Γοεμανς ανδ Δαιδ Π. Ωιλλιαμςον.
>  Ιμπροεδ αππροξιματιον αλγοριτημς φορ μαξιμυμ ςυτ ανδ ςατιςφιαβιλιτψ
> προβλεμς υςινγ ςεμιδεφινιτε προγραμμινγ.
>  Θ. Α῝Μ, 42(6):1115–1145, Νοεμβερ 1995.
>  ηττπ://ωωω-ματη.μιτ.εδυ/῀γοεμανς/ΠΑΠΕΡΣ/μαξςυτ-ϑαςμ.πδφ.
>
>
>   * [9]
> Μαρϰ Ηαλλ, Ειβε Φρανϰ, Γεοφφρεψ Ηολμες, Βερνηαρδ Πφαηρινγερ, Πετερ
> Ρευτεμανν,
> ανδ Ιαν Η. Ωιττεν.
>  Τηε ΩΕΚΑ δατα μινινγ ςοφτωαρε: Αν υπδατε.
>  ΣΙΓΚΔΔ Εξπλορ. Νεωςλ., 11(1):10–18, Νοεμβερ 2009.
>
> ηττπ://ωωω.ςιγϰδδ.οργ/ςιτες/δεφαυλτ/φιλες/ιςςυες/11-1-2009-07/π2῞11ν1.πδφ.
>
>
>   * [10]
> Θαμες Θοψςε.
>  ϒλψςςες.
>  1922\.
>  ηττπς://αρςηιε.οργ/δεταιλς/υλψςςες04300γυτ. Αςςεςςεδ Αυγυςτ
> 18, 2013.
>
>
>   * [11]
> Μοςηε Κοππελ, Ναοτ Αϰια, Ιδαν Δερςηοωιτζ, ανδ Ναςηυμ Δερςηοωιτζ.
>  ϒνςυπεριςεδ δεςομποςιτιον οφ α δοςυμεντ ιντο αυτηοριαλ ςομπονεντς.
>  Ιν Προςεεδινγς οφ τηε 49τη Αννυαλ Μεετινγ οφ τηε Αςςοςιατιον φορ
> ῝ομπυτατιοναλ Λινγυιςτιςς: Ηυμαν Λανγυαγε Τεςηνολογιες - ῞ολυμε 1, ΗΛΤ ´11,
> παγες 1356–1364, Στρουδςβυργ, ΠΑ, ϒΣΑ, 2011. Αςςοςιατιον φορ ῝ομπυτατιοναλ
> Λινγυιςτιςς.
>  ηττπ://ωωω.αςλωεβ.οργ/αντηολογψ/Π11-1136.πδφ.
>
>
>   * [12]
> Ανδρεω ΜςΔοναλδ, Σαδια Αφροζ, Αψλιν ῝αλιςϰαν, Αριελ Στολερμαν, ανδ Ραςηελ
> Γρεενςταδτ.
>  ϒςε φεωερ ινςτανςες οφ τηε λεττερ `ὶ´´: Τοωαρδ ωριτινγ ςτψλε
> ανονψμιζατιον.
>  Ιν Προςεεδινγς οφ τηε 12τη Πριαςψ Ενηανςινγ Τεςηνολογιες
> Σψμποςιυμ (ΠΕΤΣ 2012). Σπρινγερ, Θυλψ 2012.
>  ηττπς://ωωω.ςς.δρεξελ.εδυ/῀ςα499/παπερς/ανονψμουτη.πδφ.
>
>
>   * [13]
> Εφςτατηιος Σταματατος.
>  Α ςυρεψ οφ μοδερν αυτηορςηιπ αττριβυτιον μετηοδς.
>  Θ. Αμ. Σος. Ινφ. Σςι. Τεςηνολ., 60(3):538–556, Μαρςη 2009.
>
> ηττπ://ωωω.ιςςδ.αεγεαν.γρ/λεςτυρερς/ςταματατος/παπερς/ςυρεψ.πδφ.
>
>
>

>
>
>
>
>
>


# Unveiling the Nature of Infrared Bright, Optically Dark Galaxies with Early
JWST Data

L. Barrufet1, P. A. Oesch1,2, A. Weibel1, G. Brammer2, R. Bezanson3, R.
Bouwens4, Y. Fudamoto5,6<EMAIL_ADDRESS>V. Gonzalez7, G.
Illingworth8, K. E. Heintz2, B. Holden8, I. Labbe9, D. Magee8, R. P. Naidu10,
E. Nelson11 M. Stefanon12,13, R. Smit14, P. van Dokkum15, J. R. Weaver16, C.
C. Williams17
1Department of Astronomy, University of Geneva, Chemin Pegasi 51, 1290
Versoix, Switzerland
2Cosmic Dawn Center (DAWN), Niels Bohr Institute, University of Copenhagen,
Jagtvej 128, København N, DK-2200, Denmark
3Department of Physics and Astronomy and PITT PACC, University of Pittsburgh,
Pittsburgh, PA 15260, USA
4Leiden Observatory, Leiden University, PO Box 9500, 2300 RA Leiden, The
Netherlands
5Waseda Research Institute for Science and Engineering, Faculty of Science and
Engineering, Waseda University, 3-4-1 Okubo, Shinjuku, Tokyo 169-8555, Japan
6National Astronomical Observatory of Japan, 2-21-1, Osawa, Mitaka, Tokyo,
Japan
7Departamento de Astronomia, Universidad de Chile, Camino del Observatorio
1515, Las Condes, Santiago 7591245, Chile
8Department of Astronomy and Astrophysics, University of California, Santa
Cruz, CA 95064, USA
9Centre for Astrophysics and Supercomputing, Swinburne University of
Technology, Melbourne, VIC 3122, Australia
10Center for Astrophysics $|$ Harvard & Smithsonian, 60 Garden Street,
Cambridge, MA 02138, USA
11Department for Astrophysical and Planetary Science, University of Colorado,
Boulder, CO 80309, USA
12Departament d’Astronomia i Astrofisica, Universitat de Valencia, C. Dr.
Moliner 50, E-46100 Burjassot, Valencia, Spain
13Unidad Asociada CSIC "Grupo de Astrofisica Extragalactica y Cosmologia"
(Instituto de Fisica de Cantabria - Universitat de Valencia)
14Astrophysics Research Institute, Liverpool John Moores University, 146
Brownlow Hill, Liverpool L3 5RF, UK
15Astronomy Department, Yale University, 52 Hillhouse Avenue, New Haven, CT
06511, USA
16Department of Astronomy, University of Massachusetts, Amherst, MA 01003, USA
17 Steward Observatory, University of Arizona, 933 N. Cherry Avenue, Tucson,
85721, USA

(Accepted XXX. Received YYY; in original form ZZZ)

###### Abstract

Over the last few years, both ALMA and Spitzer/IRAC observations have revealed
a population of likely massive galaxies at $z>3$ that was too faint to be
detected in HST rest-frame ultraviolet imaging. However, due to the very
limited photometry for individual galaxies, the true nature of these so-called
HST-dark galaxies has remained elusive. Here, we present the first sample of
such galaxies observed with very deep, high-resolution NIRCam imaging from the
Early Release Science Program CEERS. 33 HST-dark sources are selected based on
their red colours across 1.6 $\mu$m to 4.4 $\mu$m. Their physical properties
are derived from 12-band multi-wavelength photometry, including ancillary HST
imaging. We find that these galaxies are generally heavily dust-obscured
($A_{V}\sim 2$ mag), massive ($\log(M/M_{\odot})\sim 10$), star-forming
sources at $z\sim 2-8$ with an observed surface density of $\sim 0.8$
arcmin-2. This suggests that an important fraction of massive galaxies may
have been missing from our cosmic census at $z>3$ all the way into the
reionization epoch. The HST-dark sources lie on the main-sequence of galaxies
and add an obscured star formation rate density (SFRD) of
$\mathrm{1.3^{+1.6}_{-1.0}\times 10^{-3}M_{\odot}/yr/Mpc^{3}}$ at $z\sim 6$,
similar to previous estimates. Our analysis shows the unique power of JWST to
reveal this previously missing galaxy population and to provide a complete
census of galaxies at $z=2-8$ based on rest-frame optical imaging.

###### keywords:

Galaxies: high-redshift. Infrared: galaxies

††pubyear: 2022††pagerange: Unveiling the Nature of Infrared Bright, Optically
Dark Galaxies with Early JWST Data–Unveiling the Nature of Infrared Bright,
Optically Dark Galaxies with Early JWST Data

##  1 Introduction

By $z>3$, the rest-frame optical light redshifts out of the view of the Hubble
Space Telescope (HST). Therefore, our understanding of ‘normal’ star-forming
galaxies such as Lyman Break Galaxies (LBGs) that were identified in HST
strongly relies on rest-frame ultra-violet (UV) observations. This can result
in an incomplete galaxy census at earlier times, due to UV-faint galaxy
populations such as quiescent or dust-obscured sources.

While a small number of extremely obscured, star-forming galaxies (e.g., sub-
millimeter galaxies; SMGs) have been known to exist at $z>4$ for a long time,
such galaxies are 100$\times$ less common than LBGs (with a sky density of
only 0.01 arcmin-2; e.g., Riechers et al. 2013; Marrone et al. 2018). Hence,
they only provide a limited contribution to the integrated SFRD. Similarly,
massive ($\mathrm{>10^{11}M_{\odot}}$) quiescent objects have only been
identified to $z\sim 3-4$ (e.g. Glazebrook et al., 2017; Tanaka et al., 2019;
Valentino et al., 2020; Carnall et al., 2020), but they remain relatively rare
such that degree-scale survey volumes are needed for their detection.

Nevertheless, due to our limited sensitivity in rest-frame optical light at
$z>3$, it is still possible that the early cosmic star-formation and stellar
mass densities are dominated by less extreme, dusty star-forming galaxies
(DSFGs) or quiescent galaxies, which have been missing from rest-UV datasets.
The key question is: how common are such red sources at $z>3$?

Over the last few years, evidence has been mounting that, indeed, a population
of obscured, star-forming galaxies with less extreme SFRs than SMGs could be
quite common in the early Universe. These sources have been undetected in
current HST surveys, due to their very red colors between $1.6\,\mu$m and
$4.5\,\mu$m (e.g., Franco et al. 2018; Wang et al. 2019), and they have thus
been named HIEROs (HST to IRAC extremely red objects), H-dropouts, HST-dark,
or HST-faint.

Very red dusty sources exist even below the sensitivity limits in IRAC, making
them only detectable by ALMA. However, previous studies with ALMA relied on
small areas of sky, thus reporting only one or two galaxies in some cases
(e.g., Williams et al., 2019; Fudamoto et al., 2021). Hence, HST-dark sources
are typically selected in Spitzer/IRAC imaging or as serendipitous sources in
ALMA sub-mm continuum data (e.g., Huang et al., 2011; Caputi et al., 2015;
Stefanon et al., 2015; Wang et al., 2016; Franco et al., 2018; Alcalde
Pampliega et al., 2019; Wang et al., 2019; Yamaguchi et al., 2019; Williams et
al., 2019; Umehata et al., 2020; Sun et al., 2021; Fudamoto et al., 2021;
Manning et al., 2022). However, due to the limited photometric information
available, their spectral energy distributions (SEDs) and redshifts remained
uncertain. Nevertheless, the general consensus in the literature is that the
vast majority of these galaxies are heavily dust-obscured, star-forming
sources at $z\sim 3-6$ that dominate the high-mass end of the galaxy
population at these early cosmic times (e.g., Alcalde Pampliega et al., 2019;
Wang et al., 2019).

Accurately constraining the number densities and the physical properties of
such sources is critical for our understanding of early galaxy build-up. It is
conceivable that dust-obscured, HST-undetected sources dominate the SFRD at
$z>4$ unlike what is typically assumed (e.g., Williams et al., 2019; Zavala et
al., 2021). However, the most recent 2mm number counts as well as IR
luminosity functions measured from ALMA continuum observations out to $z\sim
4-7$ indicate a less extreme scenario (see Casey et al., 2018; Gruppioni et
al., 2020, Barrufet in preparation). Nevertheless, it has become clear that
our census of galaxies at $z>3$ is indeed far from complete, leaving a
critical gap in our understanding of early galaxy build-up.

Finding and characterizing this missing galaxy population is a key goal of
extragalactic astronomy in order to inform theoretical models of galaxy
evolution. In this paper, we address this important question with some of the
first JWST NIRCam imaging from the Early Release Science (ERS) program CEERS
(Finkelstein et al., 2022). This allows us to identify HST-dark galaxies with
the same selection that has been used in the past for HST+Spitzer samples, but
now with significantly measured multi-wavelength photometry in several bands
to derive their redshifts and physical properties reliably for the first time.
This follows previous analyses with JWST data that have already revealed
massive galaxies at $z>7$ (see Santini et al., 2022; Labbe et al., 2022). In a
separate paper, we also analyze the morphologies of bright F444W selected
sources (Nelson et al. 2022, submitted).

This paper is organized as follows. In Section 2 we describe the observational
data and sample selection. Section 3 outlines our spectral energy distribution
modelling and the variety of physical properties that we find for HST-dark
galaxies. In Section 4, we present our results and compare our galaxy sample
to SMGs, before we derive their contribution to the CSFRD in Section 5 and
finish with a summary and conclusions in Section 6.

Magnitudes are given in the AB system (e.g., Oke & Gunn, 1983), and where
necessary we adopt a Planck Collaboration et al. (2016) cosmology.

##  2 Observations and Sample Selection

###  2.1 CEERS observations and merged HST-JWST catalogue.

This analysis is based on the first four NIRCam pointings from the Early
Release Science (ERS) Program CEERS (PI: Finkelstein; PID: 1345; see
Finkelstein et al., in prep.). These images include exposure sets of short-
and long-wavelength filters: F115W+F277W, F115W+F356W, F150W+F410M, and
F200W+F444W. The typical integration time for each filter set is 0.8 hr.
Hence, the F115W images obtained double the exposure time compared to the
remaining filters.

Figure 1: Colour-magnitude diagram (F444W vs F160W-F444W) used to select the
HST-dark galaxies. The dashed line shows the location of $H_{160}=27$ mag,
which corresponds to the average $3\sigma$ limit of the shallowest parts of
the HST data. The black line marks the adopted colour cut
$\mathrm{H-[4.5]>2.3}$ and the red dots correspond to the galaxies included in
our sample, HST-IRAC red sources. The black dots represent the density of
sources from the full NIR selected CEERS catalogue with detections in the
F444W band ($\mathrm{\sim 12,000}$ sources).

The data were retrieved from the MAST archive and were then WCS-matched and
combined using the publicly available
grizli111https://github.com/gbrammer/grizli/ pipeline. Most notably, the
images are aligned to the Gaia DR3 catalogs. For details on the reduction, we
refer the reader to Brammer et al., in prep. In the following, we work with
images at 40mas pixel scale.

The 5$\sigma$ depth of the new JWST data are measured in small circular
apertures of 0$\aas@@fstack{\prime\prime}$32 diameter. They range from 28.1
mag in the medium band filter F410M to 28.9 mag in the wide filters at
$2-3.5\,\mu$m (see also Table 1 in Naidu et al. 2022).

Figure 2:  Examples of four HST-dark galaxies at $\mathrm{z\sim 2.5,\ 4,\ 6,\
8}$ that represent the variety of physical properties and morphologies that we
find in our sample (see also Nelson et al. 2022). (Top panels) Postage stamps
in the optical/NIR for four examples of HST-dark galaxies, where the sources
show no significant detections in the HST-filters. The first row shows HST
images with the bands F606W, F814W, F850LP, F125W, F140W, F160W and the
Spitzer bands IRAC1, IRAC2 at $\mathrm{3.6and4.5\mu m}$ respectively. The
second row shows new NIRCam JWST images: F115W, F150W, F200W, F277W, F356W,
F410M, F444W. The extraordinary improvement in both depth and resolution of
JWST at $3-5\,\mu$m compared to the IRAC images is obvious. Most of these
sources would not have been easily selected in previous data and they would
have potentially been missing from our cosmic census. (Bottom panels)
Posterior spectral energy distributions (yellow; median with 16-84 percentile
uncertainties) with the photometry (red dots). Orange points indicate the
expected model fluxes from the posterior SEDs. The probability distribution
function of the photometric redshift is shown in the lower right part of the
panel. In general, HST-dark sources are found to be dust-obscured, massive,
star-forming galaxies at $z\sim 2-8$.

Since the CEERS survey covers the previous CANDELS/EGS field, we also include
the available ancillary data from HST. In particular, we include a re-
reduction of all the ACS and WFC3/IR imaging taken in these fields in the five
filters: F606W, F814W, F125W, F140W, and F160W. The 5$\sigma$ depths measured
in the 0$\aas@@fstack{\prime\prime}$32 circular apertures for the HST data
range from 28.4 mag in F606W to 27.0 mag in F140W. While the HST NIR data are
shallower, they still provide additional multi-wavelength constraints. Most
importantly using the F160W filter allows us to use HST-dark galaxy selections
analogous to previous HST+Spitzer searches (see next section).

The multi-wavelength photometric catalogs were derived with SExtractor (Bertin
& Arnouts, 1996). Sources were detected in the F444W image and photometry in
other bands was measured in dual-image mode. In particular, we measure the
fluxes in circular apertures of 0$\aas@@fstack{\prime\prime}$32 diameter,
before correcting them to total using the AUTO flux measurement in the F444W
image, which is based on Kron-like apertures. The input rms maps were
normalized to match the true flux uncertainties measured in random sky
positions. To account for systematic uncertainties in the photometry we apply
an error floor of 5% for all flux measurements.

The final catalog includes about 17,000 sources with photometry in 12 bands,
spanning 0.6 to 4.4 $\mu$m.

###  2.2 HST-dark Galaxy Color Selection with JWST

Previous HST-dark galaxy studies were based on analyses of HST+Spitzer/IRAC
data (e.g., Wang et al., 2016, 2019; Alcalde Pampliega et al., 2019). These
selections were limited by two factors: (1) the very broad PSF of the IRAC
data of $\sim 1.7$ ″ and (2) the limited depth of available IRAC imaging.
Hence, typical HST-dark samples were limited to very bright sources with
$\mathrm{[4.5]\lesssim 24.5mag}$. Now, with the unparalleled high-resolution,
deep JWST data, we can push such analyses to much fainter sources.
Nevertheless, in this first paper, we focus on sources that are identified
using very similar selection methods as previous HST+Spitzer samples in order
to reveal their nature.

We use the ‘traditional’ HST-dark color cut of $\mathrm{F160W-[4.5]>2.3}$ as
adopted in the literature (e.g., Caputi et al., 2012; Wang et al., 2016). This
color cut is specifically designed to identify dusty galaxies at $z\sim 3-6$.
Indeed, previous analyses have tentatively confirmed that such red HST-dark
sources lie at these redshifts on average. However, the available photometry
for these sources is very poor, with significant detections only in the IRAC
3.6 $\mu$m and 4.5 $\mu$m bands, and mostly upper limits at HST wavelengths.

We thus apply a color cut of F160W-F444W$>2.3$ to our multi-wavelength catalog
in order to select HST-dark analog sources. After also limiting the sample to
$\mathrm{F160W>27mag}$ and S/N$>$5 in the F444W filter, we identify a total of
39 galaxies in the $\sim 40$ arcmin2 down to a F444W=26.4 magnitude (see
Figure 1). We remove an additional five galaxies that do not have full
coverage in all NIRCam filters. And one more galaxy for which the photometry
is extremely noisy, likely due to remaining detector artefacts. Hence, our
final sample consists of 33 HST-dark galaxies.

##  3 SED fitting and variety of HST-dark galaxies.

We have derived the physical properties of the 33 HST-dark galaxies by using
Bayesian Analysis of Galaxies for Physical Inference and Parameter EStimation
(BAGPIPES Carnall et al., 2018). This code generates complex model galaxy
spectra to fit these to arbitrary combinations of spectroscopic and
photometric data using the MultiNest nested sampling algorithm (Feroz et al.,
2009).

The assumptions by Wang et al. (2016, 2019) have proven effective to
characterize HST-dark galaxies and this work follows the same SED modeling.
Specifically, we use a delayed star formation history (SFH) with an e-folding
time from $\mathrm{\tau=0.1-9Gyr}$ and a uniform prior. Hence, we include
relatively short bursts as well as, effectively, a constant SFH.

The stellar models are based on the 2016 updated version of the Bruzual &
Charlot (2003) library with a Kroupa initial mass function. We allow for a
range of metallicities from 0.2 to 2.5 $Z_{\odot}$. Nebular continuum and line
emission are added self-consistently based on the photoionization code CLOUDY
(Ferland et al., 1998) with the ionization parameter (log U) set as $\log
U=-2$. Finally, we adopt Calzetti et al. (2000) dust attenuation allowing for
very heavily attenuated SEDs by setting a range $\mathrm{A_{v}=0-6\,mag}$,
again using a uniform prior.

The above settings adequately reproduce the range of SEDs of our sample. We
find good fits for all our sources.

Figure 2 shows postage stamps and SED fits of four of these sources at
$\mathrm{z\sim 2.5,\ 4,\ 6,\ 8}$. They represent the variety of galaxy SEDs
that are identified as HST-dark sources as well as in terms of morphology. We
also show the cutouts of previous IRAC images from the S-CANDELS program
(Ashby et al., 2015). The outstanding increase in depth and spatial resolution
at $3-5\,\mu$m can be appreciated. With the exception of a handful of our
sample galaxies, none of these would have been easily selected based on
previously available data. Hence, the majority of our sample was completely
missing from our previous cosmic census of HST+Spitzer. In the following, we
use our sample of 33 sources to analyze the physical properties of HST-dark
galaxies. The main parameters are summarised in Table 1 in the appendix.

Figure 3: Histograms of the main physical properties analysed: photometric
redshift, SFR (first row), stellar masses and dust attenuation (second row).
The histograms show the generally dusty high-z nature of the sample with most
of the galaxies at $z>2.5$ and a median attenuation of $\mathrm{A_{v}\sim 2}$.
The SFRs are moderate, mostly at $\mathrm{SFR<100M_{\odot}/year}$ whereas the
stellar masses are large with $\mathrm{log(M_{*}/M_{\odot})\sim 10}$. Our
analysis thus shows that previous stellar mass function estimates at $z>3$
were likely underestimated at high masses.

##  4 HST-dark galaxies: The high-z extension of SMGs?

###  4.1 HST-dark galaxy properties

Thanks to the very sensitive NIRCam data, we can now – for the first time –
measure accurate photometric redshifts and SEDs for HST-dark galaxies. Hence,
we can now derive significantly more reliable physical parameters such as SFRs
and stellar masses compared to previous IRAC-selected HST-dark sources. In
this section, we thus present the physical properties derived from SED fitting
for our final sample of 33 sources (see Figures 2 and 3).

We find that an HST-dark galaxy selection includes a variety of sources with
the vast majority of the sample being high redshift dust-obscured, massive,
star-forming systems. The photometric redshifts are accurate for the galaxies
at $z>2.5$ with small uncertainties ($\Delta z\sim 0.1-0.2$) with only one
solution in the redshift posterior distribution function (see table 1 and
Figure 2). The subset of galaxies that is placed at $z<2.5$ presents more
considerable uncertainties and typically shows two possible redshift
solutions. As expected, we confirm the dusty nature of these colour-selected
HST-dark galaxies. We have allowed the dust attenuation to vary broadly in the
modelling, but the values only reach up to $\mathrm{A_{v}=4mag}$ for almost
all HST-dark galaxies with a median value for the sample of $\mathrm{A_{v}\sim
2.3mag}$.

Interestingly, the sample reaches redshifts $\mathrm{z\sim 8}$ for several
galaxies. In total, we identify eleven galaxies at $\mathrm{z>6}$, i.e., in
the epoch of reionization, whereas previous HST-dark galaxies selected with
Spitzer were typically limited to $z<6$. Only one of our galaxies is in common
with the study of Labbe et al. (2022), who specifically selected the most
massive galaxies at high-redshifts in the early JWST data. However, these
sources at $z>6$ are certainly the most surprising due to their high mass and
high dust extinction, identified in only 40 arcmin2 of JWST imaging. The high
attenuation for the galaxies above $\mathrm{z>6}$ suggests the presence of
dust in the reionization epoch (see also Fudamoto et al., 2021). We refer to a
future paper for more detailed analyses of these sources.

We have also analysed both the stellar masses and the SFRs from the SED
fitting. The SFRs are moderate, when compared to SMGs, with $\mathrm{97\%}$ of
our sample presenting $\mathrm{SFR<100M_{\odot}/yr}$. The stellar masses of
our sample are $\mathrm{log(M_{*}/M_{\odot})>9}$ and for the eleven galaxies
at $\mathrm{z>6}$ we find $\mathrm{log(M_{*}/M_{\odot})=9.8\pm 0.3}$.

Figure 4:  HST-dark galaxies compared with submillimetre galaxies (SMGs). (Top
panel) The redshift histogram shows that HST-dark galaxies (red line) lie
mostly at $\mathrm{z>2.5}$ and extend to higher redshift compared to SMGs
(orange bars). (Bottom panel) Photometric redshifts against star formation
rates. On average, HST-dark galaxies lie at higher redshift and have
significantly lower SFRs (orange squares; ALESS sample from da Cunha et al.
(2015)). The significantly smaller error bars in the photometric redshifts are
due to the improvement coming from the deep, high-resolution JWST/NIRCam data.

###  4.2 HST-dark galaxies compared to SMGs.

Given their dusty nature, it is interesting to compare the HST-dark galaxies
with a classical SMGs sample, in order to put them in a broader context of
massive high-z dusty galaxies. We find that HST-dark galaxies have higher
redshift compared to SMGs (see the histogram in Figure 4). Interestingly and
probably not unexpectedly, HST-dark galaxies have lower SFRs than SMGs.
Similarly, we also find HST-dark galaxies have more than an order of magnitude
lower stellar masses than the massive SMGs that weigh $\mathrm{M_{*}\sim
10^{12}M_{\odot}}$.

In Figure 5, we further show the location of HST-dark and SMGs with respect to
the main sequence (MS) of galaxies at $z\sim 4$. As can be seen, HST-dark
sources at these redshifts follow the MS. They populate its lower mass end
compared to both SMGs and IRAC-selected HST-Dark galaxies. Specifically, we
have compared our results with two previous studies of dusty galaxies: The
ALESS sample from da Cunha et al. (2015) shows a larger number of starburst
galaxies that lie significantly above the MS. In between our sample and SMGs,
we find more massive HST-dark galaxies from Wang et al. (2019). This
difference could be due to the selection criteria and lower mass HST-dark
galaxies could have been missed by using Spitzer selection.

In a nutshell, HST-dark galaxies appear to be the higher redshift, lower SFR
extension of SMGs. They have likely just not been detected in previous sub-
millimeter single-dish studies due to their faintness. The synergy between
JWST and ALMA data will be needed to confirm this rather interesting outcome.

Figure 5: Star formation rate against stellar mass for our sample
corresponding to $3<z<5$ (red diamonds) and the total sample (empty diamonds).
The black line shows the MS at $z=4$ from Schreiber et al. (2015). HST-dark
galaxies lie on the MS as other optically dark massive galaxies selected with
Spitzer data from Wang et al. (2019) (yellow dots). SMGs from da Cunha et al.
(2015) are above the MS (orange squares).

##  5 Cosmic Star Formation Rate: is HST missing a significant contribution at
high-z?

The cosmic star formation history is key to explaining the mass assembly and
the metal enrichment of the galaxy population both vital to understanding
galaxy evolution. Although the Cosmic Star Formation Rate Density (CSFRD) is
well constrained up to $z\sim 8$ for UV-selected galaxies (i.e., LBGs), the
contribution from obscured, IR rest-frame selected galaxies is still highly
uncertain and only poorly constrained at $z>3$ (e.g., Casey et al., 2014,
2019; Zavala et al., 2021). A central question is: when does the dust-obscured
star formation start to dominate the CSFRD? And, more specifically, what is
the contribution to the CSFRD from HST-dark galaxies? In this Section, we
approach this question and put our results in the context of the most relevant
high-$z$ galaxy populations.

We calculate the SFRD for our sample of HST-dark galaxies at $z>3.5$ in three
redshift bins with the same number of galaxies centered at $z=3.7,5.7,7.2$
respectively. We find a contribution of $\mathrm{SFRD=1.3^{+2.5}_{-0.5}\times
10^{-3}M_{\odot}/yr/Mpc^{3}}$ at $z\sim 4$ that remains constant within the
errorbars up to $z\sim 7$. Specifically, we measure $\mathrm{SFRD\sim
1.3^{+1.6}_{-1.0}\times 10^{-3}M_{\odot}/yr/Mpc^{3}}$ at $z\sim 6$, and
$\mathrm{SFRD\sim 3.2^{+4.1}_{-2.6}\times 10^{-3}M_{\odot}/yr/Mpc^{3}}$ by
$z\sim 7$.

Figure 6 shows that our results have lower SFRD than previously selected
H-dropout galaxy samples at $z\sim 4-5$, but are still in agreement with the
uncertainty. There are no current observations of HST-dark galaxies at $z>6$.
Our results however suggest that the SFRD of HST-dark galaxies remains
constant or even larger SFRD at $z\sim 7$. Our redshift bins currently only
contain six sources each. Therefore, future studies with larger samples are
needed to corroborate this result.

Our analysis shows an HST-dark contribution to the SFRD similar as SMGs at
$4<z<5$, but larger at $z>5$. An evident explanation is that the number of
SMGs drops at $z>5$ (da Cunha et al., 2015). Even if their SFRs are an order
of magnitude higher than the average of HST-dark galaxies, their SFRD is lower
at high-$z$.

In order to detect fainter dusty star forming galaxies, DSFGs, at higher
redshifts than SMGs, Casey et al. (2021) have recently conducted the ALMA MORA
survey using 2mm continuum observations over 184 $\mathrm{arcmin^{2}}$. Based
on these data, they identify galaxies with lower SFRs than the extreme SMGs,
but with a higher number density at $\mathrm{z\sim 4}$, probably due to a
selection effect. Our study shows that HST-dark galaxies contribute an order
of magnitude less than 2mm-selected DSFGs. The reason probably is the lower
SFRs and this suggests that fainter dusty galaxies could have been missed by
shallow ALMA surveys at $z\sim 4-7$.

Our study, and previous studies mentioned above, indicate that the dusty
population at high-$z$ is larger than what was expected. In particular, more
extensive samples are needed from deeper surveys and multi-wavelength
observations to understand how the disparate objects seen today relate and
contribute to the overall cosmic SFRD. A very important result from our
analysis is that HST-dark galaxies significantly contribute to the CSFRD, in
particular at $z>5$.

Figure 6:  Star formation rate density against redshift for the HST-dark
galaxies (red diamonds) compared to a compilation of previous studies. The
orange area represents the dust-obscured star formation history from Zavala et
al. (2021). The total SFRD as a function of redshift is shown as the black
dashed line from Madau & Dickinson (2014). The contribution of our HST-dark
galaxies remains almost constant from redshift $3.5<z<6$. The HST-dark
contribution is similar to SMGs at $z\sim 5-6$ (orange squares; da Cunha et
al. 2015). Our SFRD is in agreement with previous HST-dark studies studies at
$z\sim 4-6$ (orange dots; Wang et al. 2019) but we extend the result to $z\sim
7$. DSFGs from Casey et al. 2021 show similar SFRD at $z\sim 6$ as our JWST
selected HST-dark galaxies.

##  6 Summary and Conclusions

In this work, we have exploited some of the first JWST/NIRCam data to identify
red high-$z$ galaxies that were missed by HST’s rest-frame UV selections.
These HST-dark galaxies were selected based on $H-F444W$ colors. We studied 33
galaxies in the four ERS CEERS NIRCam pointings, and used SED fitting to
derive their properties and put them into context by comparison to other dusty
high-redshift galaxy populations.

The main results are summarized as follows:

  * •

The average properties of HST-dark galaxies show their high-$z$ dusty nature
with 82 % of the sample at $2\lesssim z\lesssim 8$. Surprisingly, we identify
eleven galaxies at $z>6$ with a median dust attenuation of $\mathrm{A_{v}\sim
1.8\pm 0.4mag}$. This indicates that massive, heavily obscured galaxies have
been missing from our cosmic census in the epoch of reionization.

  * •

We find that HST-dark galaxies have moderate star formation rates compared to
SMGs with $\mathrm{SFR<100,M_{\odot}/yr}$, but they extend to higher
redshifts. HST-dark galaxies lie on the MS of galaxies with lower masses than
SMGs, but still with a median of $\mathrm{log_{10}(M_{*}/M_{\odot})=9.8\pm
0.3}$ for the eleven galaxies at $z>6$.

  * •

We find that the contribution of HST-dark galaxies to the cosmic star
formation rate density (CSFRD) is $1.3^{+1.6}_{-1.0}\times
10^{-3}M_{\odot}/\mathrm{yr}/\mathrm{Mpc}^{3}$ at $z=5.7$ remaining almost
constant over redshift $3.5<z<7$. Our SFRD results are in broad agreement with
previous HST-dark galaxy studies, but the JWST-selected sample indicates that
HST-dark galaxies existed up to higher redshifts than previously assumed, i.e.
$z\sim 8$.

Thanks to the outstanding improvement of NIRCam data over the previous
Spitzer/IRAC imaging, this work provides more reliable photometric redshifts,
SFRs, and stellar masses compared to previous HST-dark galaxy studies. Hence,
we can finally conclude that HST-dark galaxies have a high-z nature and that
they are numerous at $z>2.5$, with a surface density of $\sim 0.8$ arcmin-2.
The high number of dusty galaxies with $z>6$ suggest that the presence of dust
is still relevant in the reionization epoch.

The SFRD contribution of our sample of HST-dark galaxies is in agreement with
previous H-dropout studies at $z<6$. We also find a significant SFRD at $z\sim
7$ that needs to be corroborated with larger samples of HST-dark galaxies. But
already now, with these first NIRCam images, JWST has demonstrated its power
to identify more complete samples of rest-frame optically selected galaxies
out to very high redshifts. We have shown that dusty galaxies at $z>3$ are
likely to be far more numerous than expected and that such HST-dark galaxies
provide an important, but hitherto underestimated contribution to the high
mass end of the galaxy mass function and the cosmic SFRD up to $z\sim 7-8$.

## Acknowledgements

We acknowledge support from: the Swiss National Science Foundation through
project grant 200020_207349 (LB, PAO, and AW). The Cosmic Dawn Center (DAWN)
is funded by the Danish National Research Foundation under grant No. 140. YF
acknowledge support from NAOJ ALMA Scientific Research Grant number 2020-16B.
VG acknowledges support by the ANID BASAL projects ACE210002 and FB210003 KEH
acknowledges support from the Carlsberg Foundation Reintegration Fellowship
Grant CF21-0103. RS acknowledges a STFC Ernest Rutherford Fellowship
(ST/S004831/1) MS acknowledges support from the CIDEGENT/2021/059 grant, and
from project PID2019-109592GB-I00/AEI/10.13039/501100011033 from the Spanish
Ministerio de Ciencia e Innovacion - Agencia Estatal de Investigacion. RJB and
MS acknowledge support from NWO grant TOP1.16.057.

Cloud-based data processing and file storage for this work is provided by the
AWS Cloud Credits for Research program.

This work is based on observations made with the NASA/ESA/CSA James Webb Space
Telescope. The data were obtained from the Mikulski Archive for Space
Telescopes at the Space Telescope Science Institute, which is operated by the
Association of Universities for Research in Astronomy, Inc., under NASA
contract NAS 5-03127 for JWST. These observations are associated with program
# 1345.

Facilities: JWST, HST

Software: matplotlib (Hunter, 2007), numpy (Oliphant, 2015), scipy (Virtanen
et al., 2020), jupyter (Kluyver et al., 2016), Astropy (Astropy Collaboration
et al., 2013, 2018), grizli (Brammer, 2018), SExtractor (Bertin & Arnouts,
1996), bagpipes (Carnall et al., 2018)

## Data Availability

The calibrated (stage 2) data used here are available Mikulski Archive for
Space Telescopes (https://mast.stsci.edu). Further data products are available
from the authors upon reasonable request.

## References

  * Alcalde Pampliega et al. (2019) Alcalde Pampliega B., et al., 2019, ApJ, 876, 135
  * Ashby et al. (2015) Ashby M. L. N., et al., 2015, ApJS, 218, 33
  * Astropy Collaboration et al. (2013) Astropy Collaboration et al., 2013, A&A, 558, A33
  * Astropy Collaboration et al. (2018) Astropy Collaboration et al., 2018, AJ, 156, 123
  * Bertin & Arnouts (1996) Bertin E., Arnouts S., 1996, A&AS, 117, 393
  * Brammer (2018) Brammer G., 2018, Gbrammer/Grizli: Preliminary Release, Zenodo, doi:10.5281/zenodo.1146905
  * Bruzual & Charlot (2003) Bruzual G., Charlot S., 2003, MNRAS, 344, 1000
  * Calzetti et al. (2000) Calzetti D., Armus L., Bohlin R. C., Kinney A. L., Koornneef J., Storchi-Bergmann T., 2000, ApJ, 533, 682
  * Caputi et al. (2012) Caputi K. I., et al., 2012, ApJ, 750, L20
  * Caputi et al. (2015) Caputi K. I., et al., 2015, ApJ, 810, 73
  * Carnall et al. (2018) Carnall A. C., McLure R. J., Dunlop J. S., Davé R., 2018, MNRAS, 480, 4379
  * Carnall et al. (2020) Carnall A. C., et al., 2020, MNRAS, 496, 695
  * Casey et al. (2014) Casey C. M., Narayanan D., Cooray A., 2014, Phys. Rep., 541, 45
  * Casey et al. (2018) Casey C. M., et al., 2018, ApJ, 862, 77
  * Casey et al. (2019) Casey C., et al., 2019, BAAS, 51, 212
  * Casey et al. (2021) Casey C. M., et al., 2021, ApJ, 923, 215
  * Ferland et al. (1998) Ferland G. J., Korista K. T., Verner D. A., Ferguson J. W., Kingdon J. B., Verner E. M., 1998, PASP, 110, 761
  * Feroz et al. (2009) Feroz F., Hobson M. P., Bridges M., 2009, MNRAS, 398, 1601
  * Finkelstein et al. (2022) Finkelstein S. L., et al., 2022, arXiv e-prints, p. arXiv:2207.12474
  * Franco et al. (2018) Franco M., et al., 2018, A&A, 620, A152
  * Fudamoto et al. (2021) Fudamoto Y., et al., 2021, Nature, 597, 489
  * Glazebrook et al. (2017) Glazebrook K., et al., 2017, Nature, 544, 71
  * Gruppioni et al. (2020) Gruppioni C., et al., 2020, A&A, 643, A8
  * Huang et al. (2011) Huang J. S., Zheng X. Z., Rigopoulou D., Magdis G., Fazio G. G., Wang T., 2011, ApJ, 742, L13
  * Hunter (2007) Hunter J. D., 2007, Computing In Science & Engineering, 9, 90
  * Kluyver et al. (2016) Kluyver T., et al., 2016, in Loizides F., Schmidt B., eds, Positioning and Power in Academic Publishing: Players, Agents and Agendas. pp 87 – 90
  * Labbe et al. (2022) Labbe I., et al., 2022, arXiv e-prints, p. arXiv:2207.12446
  * Madau & Dickinson (2014) Madau P., Dickinson M., 2014, ARA&A, 52, 415
  * Manning et al. (2022) Manning S. M., et al., 2022, ApJ, 925, 23
  * Marrone et al. (2018) Marrone D. P., et al., 2018, Nature, 553, 51
  * Naidu et al. (2022) Naidu R. P., et al., 2022, arXiv e-prints, p. arXiv:2207.09434
  * Oke & Gunn (1983) Oke J. B., Gunn J. E., 1983, ApJ, 266, 713
  * Oliphant (2015) Oliphant T. E., 2015, Guide to NumPy. Continuum Press
  * Planck Collaboration et al. (2016) Planck Collaboration et al., 2016, A&A, 594, A13
  * Riechers et al. (2013) Riechers D. A., et al., 2013, Nature, 496, 329
  * Santini et al. (2022) Santini P., et al., 2022, arXiv e-prints, p. arXiv:2207.11379
  * Schreiber et al. (2015) Schreiber C., et al., 2015, A&A, 575, A74
  * Stefanon et al. (2015) Stefanon M., et al., 2015, ApJ, 803, 11
  * Sun et al. (2021) Sun F., et al., 2021, arXiv e-prints, p. arXiv:2109.01751
  * Tanaka et al. (2019) Tanaka M., et al., 2019, ApJ, 885, L34
  * Umehata et al. (2020) Umehata H., et al., 2020, A&A, 640, L8
  * Valentino et al. (2020) Valentino F., et al., 2020, ApJ, 889, 93
  * Virtanen et al. (2020) Virtanen P., et al., 2020, Nature Methods, 17, 261
  * Wang et al. (2016) Wang T., et al., 2016, ApJ, 816, 84
  * Wang et al. (2019) Wang T., et al., 2019, Nature, 572, 211
  * Williams et al. (2019) Williams C. C., et al., 2019, ApJ, 884, 154
  * Yamaguchi et al. (2019) Yamaguchi Y., et al., 2019, ApJ, 878, 73
  * Zavala et al. (2021) Zavala J. A., et al., 2021, ApJ, 909, 165
  * da Cunha et al. (2015) da Cunha E., et al., 2015, ApJ, 806, 110

Table 1: HST-dark Galaxy Sample and Derived Properties ID | RA | DEC | $z$ | log Mass | log SFR | $A_{V}$
---|---|---|---|---|---|---
| [deg] | [deg] |  | [$M_{\odot}$] | [$M_{\odot}$ yr-1] | [mag]
1590† | 214.8405342 | 52.8179497 | 8.19${}^{+0.33}_{-0.23}$ | 9.91${}^{+0.16}_{-0.16}$ | 1.87${}^{+0.15}_{-0.14}$ | 2.46${}^{+0.25}_{-0.24}$
5740 | 214.9983863 | 53.0046347 | 7.96${}^{+0.09}_{-0.06}$ | 10.02${}^{+0.08}_{-0.07}$ | 2.07${}^{+0.06}_{-0.06}$ | 2.14${}^{+0.09}_{-0.10}$
4682 | 214.9568193 | 52.9731632 | 7.94${}^{+0.09}_{-0.06}$ | 9.32${}^{+0.10}_{-0.07}$ | 1.37${}^{+0.06}_{-0.06}$ | 1.36${}^{+0.10}_{-0.11}$
1942 | 214.9257399 | 52.9456743 | 7.20${}^{+0.13}_{-0.11}$ | 9.60${}^{+0.12}_{-0.10}$ | 1.64${}^{+0.08}_{-0.09}$ | 2.22${}^{+0.15}_{-0.14}$
3383 | 214.9724243 | 52.9622051 | 7.19${}^{+0.17}_{-0.10}$ | 9.98${}^{+0.13}_{-0.11}$ | 1.82${}^{+0.15}_{-0.15}$ | 2.28${}^{+0.24}_{-0.25}$
250 | 214.9233523 | 52.9255999 | 7.17${}^{+0.10}_{-0.07}$ | 9.72${}^{+0.08}_{-0.05}$ | 1.80${}^{+0.06}_{-0.04}$ | 1.84${}^{+0.07}_{-0.06}$
4309 | 214.7943808 | 52.8374062 | 7.09${}^{+0.11}_{-0.06}$ | 9.96${}^{+0.12}_{-0.13}$ | 1.79${}^{+0.06}_{-0.06}$ | 1.35${}^{+0.08}_{-0.09}$
846 | 214.8493702 | 52.8118422 | 6.66${}^{+0.08}_{-0.21}$ | 10.23${}^{+0.09}_{-0.08}$ | 1.83${}^{+0.11}_{-0.11}$ | 1.86${}^{+0.13}_{-0.11}$
967 | 214.9272223 | 52.9339071 | 6.24${}^{+0.17}_{-0.20}$ | 9.95${}^{+0.08}_{-0.11}$ | 1.26${}^{+0.16}_{-0.19}$ | 0.97${}^{+0.19}_{-0.25}$
4954 | 215.0084751 | 52.9779849 | 6.01${}^{+0.32}_{-0.18}$ | 9.21${}^{+0.12}_{-0.09}$ | 1.25${}^{+0.09}_{-0.09}$ | 1.49${}^{+0.12}_{-0.16}$
11561 | 214.8383827 | 52.8852022 | 5.99${}^{+0.26}_{-0.18}$ | 9.65${}^{+0.07}_{-0.06}$ | 1.73${}^{+0.06}_{-0.05}$ | 1.84${}^{+0.06}_{-0.07}$
6731 | 214.8867805 | 52.8553889 | 5.44${}^{+0.05}_{-0.06}$ | 10.23${}^{+0.11}_{-0.16}$ | 1.96${}^{+0.09}_{-0.10}$ | 2.47${}^{+0.15}_{-0.13}$
8173 | 214.8644764 | 52.8709989 | 5.27${}^{+0.06}_{-0.23}$ | 10.10${}^{+0.14}_{-0.13}$ | 1.95${}^{+0.10}_{-0.09}$ | 2.17${}^{+0.14}_{-0.15}$
7275 | 214.8400204 | 52.8606629 | 4.49${}^{+0.14}_{-0.12}$ | 9.56${}^{+0.06}_{-0.10}$ | 1.04${}^{+0.16}_{-0.11}$ | 1.26${}^{+0.18}_{-0.12}$
916 | 214.8918494 | 52.9338890 | 4.27${}^{+1.57}_{-2.43}$ | 10.36${}^{+0.33}_{-0.79}$ | 1.95${}^{+0.40}_{-1.20}$ | 2.85${}^{+1.33}_{-0.38}$
3471 | 214.7665658 | 52.8315351 | 4.14${}^{+0.14}_{-0.13}$ | 9.54${}^{+0.06}_{-0.07}$ | 0.98${}^{+0.12}_{-0.11}$ | 1.68${}^{+0.15}_{-0.13}$
5554 | 214.8871069 | 52.8453846 | 4.04${}^{+0.14}_{-0.34}$ | 9.50${}^{+0.05}_{-0.07}$ | 0.80${}^{+0.10}_{-0.11}$ | 1.45${}^{+0.14}_{-0.13}$
975 | 214.8550691 | 52.8130538 | 3.65${}^{+0.13}_{-0.09}$ | 9.89${}^{+0.04}_{-0.04}$ | 0.30${}^{+0.35}_{-0.22}$ | 0.52${}^{+0.43}_{-0.27}$
4207 | 214.9508216 | 52.9668752 | 3.65${}^{+0.27}_{-0.18}$ | 9.61${}^{+0.10}_{-0.13}$ | 1.00${}^{+0.16}_{-0.16}$ | 2.33${}^{+0.18}_{-0.17}$
3550 | 214.8081566 | 52.8322243 | 3.62${}^{+0.12}_{-0.07}$ | 10.17${}^{+0.04}_{-0.03}$ | 0.71${}^{+0.30}_{-0.27}$ | 1.49${}^{+0.31}_{-0.24}$
10249 | 214.8350965 | 52.8951403 | 3.50${}^{+0.03}_{-0.02}$ | 9.57${}^{+0.15}_{-0.13}$ | 1.45${}^{+0.04}_{-0.14}$ | 2.33${}^{+0.13}_{-0.18}$
2019 | 214.8296410 | 52.8207869 | 3.44${}^{+0.17}_{-0.20}$ | 10.45${}^{+0.09}_{-0.12}$ | 1.86${}^{+0.14}_{-0.12}$ | 3.03${}^{+0.18}_{-0.18}$
10517 | 214.8890479 | 52.8926261 | 3.32${}^{+0.12}_{-0.12}$ | 10.57${}^{+0.08}_{-0.13}$ | 1.96${}^{+0.13}_{-0.11}$ | 3.27${}^{+0.23}_{-0.27}$
10688 | 214.8894141 | 52.8922178 | 3.27${}^{+0.12}_{-0.08}$ | 8.97${}^{+0.09}_{-0.08}$ | 0.68${}^{+0.06}_{-0.11}$ | 0.96${}^{+0.09}_{-0.14}$
2604 | 214.9257588 | 52.9544549 | 2.56${}^{+0.13}_{-0.17}$ | 9.65${}^{+0.12}_{-0.16}$ | 0.93${}^{+0.20}_{-0.17}$ | 2.99${}^{+0.31}_{-0.30}$
204 | 214.8402522 | 52.8011218 | 2.53${}^{+0.07}_{-0.09}$ | 9.57${}^{+0.17}_{-0.17}$ | 1.23${}^{+0.11}_{-0.19}$ | 2.74${}^{+0.18}_{-0.24}$
6510 | 214.8475603 | 52.8533817 | 2.34${}^{+2.39}_{-0.56}$ | 9.55${}^{+0.87}_{-0.28}$ | 0.80${}^{+1.23}_{-0.41}$ | 4.71${}^{+0.85}_{-1.31}$
3740 | 214.7598105 | 52.8334256 | 2.06${}^{+1.71}_{-0.44}$ | 9.01${}^{+0.46}_{-0.33}$ | 0.22${}^{+0.71}_{-0.46}$ | 2.46${}^{+0.66}_{-0.44}$
6611 | 214.8559459 | 52.8546735 | 1.74${}^{+0.13}_{-0.13}$ | 9.90${}^{+0.14}_{-0.17}$ | 1.02${}^{+0.20}_{-0.16}$ | 4.08${}^{+0.37}_{-0.36}$
1187 | 214.9091023 | 52.9372295 | 1.73${}^{+0.85}_{-0.21}$ | 9.29${}^{+0.34}_{-0.23}$ | 0.42${}^{+0.47}_{-0.28}$ | 4.71${}^{+0.75}_{-1.06}$
5609 | 214.8706440 | 52.8461150 | 1.72${}^{+1.88}_{-0.28}$ | 9.64${}^{+0.54}_{-0.21}$ | 0.74${}^{+0.85}_{-0.29}$ | 3.27${}^{+0.51}_{-0.94}$
5705 | 215.0045402 | 52.9835368 | 1.63${}^{+0.12}_{-0.26}$ | 9.88${}^{+0.16}_{-0.18}$ | 1.08${}^{+0.27}_{-0.27}$ | 5.06${}^{+0.38}_{-0.35}$
9582 | 214.9315526 | 52.9210161 | 1.13${}^{+3.32}_{-0.07}$ | 9.39${}^{+0.86}_{-0.25}$ | -0.18${}^{+2.14}_{-2.16}$ | 5.46${}^{+0.39}_{-2.23}$
$\dagger$ source 2859 in Labbe et al. 2022