# The Projection-Enhancement Network

Christopher Z. Eddy
Department of Physics
Oregon State University
Corvallis, OR 97331
Austin Naylor
Department of Physics
Oregon State University
Corvallis, OR 97331
Bo Sun
Department of Physics
Oregon State University
Corvallis, OR 97331
Primary correspondence email<EMAIL_ADDRESS>
###### Abstract

Contemporary approaches to instance segmentation in cell science use 2D or 3D
convolutional networks depending on the experiment and data structures.
However, limitations in microscopy systems or efforts to prevent phototoxicity
commonly require recording sub-optimally sampled data regimes that greatly
reduces the utility of such 3D data, especially in crowded environments with
significant axial overlap between objects. In such regimes, 2D segmentations
are both more reliable for cell morphology and easier to annotate. In this
work, we propose the Projection Enhancement Network (PEN), a novel
convolutional module which processes the sub-sampled 3D data and produces a 2D
RGB semantic compression, and is trained in conjunction with an instance
segmentation network of choice to produce 2D segmentations. Our approach
combines augmentation to increase cell density using a low-density cell image
dataset to train PEN, and curated datasets to evaluate PEN. We show that with
PEN, the learned semantic representation in CellPose encodes depth and greatly
improves segmentation performance in comparison to maximum intensity
projection images as input, but does not similarly aid segmentation in region-
based networks like Mask-RCNN. Finally, we dissect the segmentation strength
against cell density of PEN with CellPose on disseminated cells from side-by-
side spheroids. We present PEN as a data-driven solution to form compressed
representations of 3D data that improve 2D segmentations from instance
segmentation networks.

_Keywords_ First keyword $\cdot$ Second keyword $\cdot$ More

##  1 Introduction

Automated computational methods are crucial for high-throughput analyses of
microscopy images, where structures of interest are tagged through staining,
endogenous expression of fluorophores, or identified through contrast methods.
The subsequent image processing, however, often requires expensive expert-
level identification [22, 21]. In the domain of cell-science, instance
segmentation, or the pixel-wise identification of each unique occurrence of an
object in an image, is essential to capture vital morphological and biological
insights, and has led to a deeper understanding of cell-heterogeneity [32],
the spatial-organization of sub-cellular components [13, 7], and phenotype
transitions in cancer [35, 10], to name a few. Deep Neural Networks (DNN) and
computer vision methods have been instramental to accomplish these tasks.

Many biomedical image analyses utilize convolutional neural networks for the
identification of objects of interest in their images, in part due to their
ability to learn and extract important features in the local receptive fields
of stacked convolutions [23, 3]. Many of such applications take advantage of
two particular architectures, including region-based networks, which propose
object regions in an image for downstream segmentation, and U-Net based
architectures, which contain an encoder-decoder style network that extracts
features and spatial information to construct object segmentations [18].

While many imaging modalities are able to acquire 3D spatial data, several
challenges exist in fully-realizing its utility. First, researchers are often
limited in 3D resolution due to toxicity or bleaching effects during imaging.
To address the issue, computational algorithms have been proposed to infer a
high-resolution 3D image from a sub-optimally sampled 3D image stack. The
traditional method utilizes deconvolution of the spatially anisotropic point-
spread-function with interpolation to overcome the insufficient axial
resolution, at the expense of errors in the deconvolution method and
additional parameters to hand-tune [8, 11]. More recently, state-of-the-art
resolution enhancing deep learning techniques have been proposed and proven
highly effective for both medical [6, 31] and microscopy data [36, 38, 34].

When high resolution 3D data is available, it often demands significant
overhead in computational time and memory requirements for instance
segmentation. Therefore, the majority of current methods do not use an end-to-
end approach on 3D data, and instead charge the deep learning networks to only
perform semantic segmentation, pixel-wise classification on 2D image slices,
and later processed downstream by seeded watershed [12, 18] or other
traditional segmentation techniques [37, 33]. When axial resolution is high
enough, a different strategy may be to use 2D instance segmentation networks
to label a 3D image using all available 2D slices [27].

Finally, in training of cell-based DNNs, few public sources of annotated
datasets for 3D imaging modalities are available in part due to the tedious
nature of annotating such data slice-by-slice. While promising semi-supervised
methods have been considered to cut the necessary manual labor costs of
annotating data [5], they may introduce unintended bias [4].

Due to all these constraints, it is desirable to achieve accurate cell
segmentation on 3D image stacks that are sparsely sampled along the axial (z)
dimension. The task is particularly challenging at high cell densities. 2D
instance segmentation networks have far less parameters involved and offer an
end-to-end solution to acquire 2D segmentations of objects. Moreover, there
are an abundance of large, readily available labeled 2D cell images through
Cell Image Library, Image Data Resource, and Kaggle which can be used to
easily train 2D networks. At the single cell level, 2D images also encode most
of the morphological quantities that provide accurate phenotype classification
[10].

Given the advantages of utilizing 2D images, it is imperative to recognize the
limitations of simple dimensional-reduction approaches. The widely used
maximum intensity projection (MIP), for instance, does not provide depth
features in order to maintain contextual information, and therefore, spatial
context is lost. In segmentation tasks, such as that shown in Figure 1B, MIP
introduces spurious overlapping objects that are occluded and results in
under-counting and poor segmentation. Other forms of projections, including
standard deviation, sum, and mean projections can each introduce their own
artifacts into the compressed representation. Instead, many researchers have
used color and depth image pairs to overcome the loss of 3D spatial cues to
perform 2D instance segmentation [26, 15, 14].

In order to assist cell segmentation in 3D images that are sparsely sampled
along the axial dimension, we develop the projection-enhancement network
(PEN). PEN is a fully convolutional module which acts as a data-driven unit to
encode spatial information of 3D microscopy images into compressed 2D
representations. As a module, PEN is placed in front of the 2D instance
segmentation network of choice and is trained concurrently to maximize the
learning objectives of the instance segmentation network. We show that in
contrast to MIP methods, PEN results in significant gains in detection and
segmentations in high-density cells in 3D cultures. We show that functionally,
PEN learns to encode depth, or becomes a low-high pass filter depending on the
training setup. We highlight the segmentation ability of PEN in cancer cells
disseminating from spheroids. Considering these results, we present PEN as an
effective tool to decrease critical computation time and provide a method to
spatially resolve 3D distributed objects in microscopy images for downstream
analyses.

##  2 Results

In order to take advantage of 2D image segmentation techniques on 3D image
stacks that have low axial resolution, we propose a data-driven model to
optimally reduce a gray scale image stack to a 2D RGB image. Our model is
inspired by the Inception module [29]. A requisite of the model design
included forming a shallow network to limit the overhead in terms of memory
and computation time, as PEN is built in line with 2D segmentation networks as
shown in Figure 1C. Cells distributed in 3D may take any orientation and vary
in shape and spatial distribution. These challenges lead us to select an
architecture of a wide network, which performs independent operations at
multiple scales that are concatenated at the output step.

Specifically, PEN consists of 3D convolutions distributed in separate
branches, as shown in Figure 1A. In each branch, a single convolutional kernel
of size K is applied to the 3D image of axial size Z without padding in the
axial dimension, and forms 3 feature maps. Following all convolutions, the
feature maps undergo ReLU activation, then batch normalization. A subsequent
convolution with kernel size of (1, 1, $Z-K$) is applied to pool the axial
features. The axial dimension is then squeezed out, and the semantic image in
each branch becomes a 2D image with RGB channels. The outputs of the branches
are then stacked, and a final 3D convolution is applied with kernel sizes of
(1, 1, $N_{branches}$) and 3 output channels. The convolution acts to pool
each branch image separately into each output color channel, followed by non-
linear ReLU activation and batch-normalization. The third spatial dimension is
then squeezed out, leaving a 2D RGB-color image that is rescaled and
normalized to be fed to the 2D segmentation network of choice.

Figure 1: The Projection-Enhancement Network (PEN). (A) Architecture of PEN to
encode 3D axial data into a 2D output image; A Z-stack 3D image is passed as
input to PEN, which is operated on by 5 different scales of conv blocks. The
outputs of each branch are 3 x H x W, which are then stacked together, and
operated on by a final conv block to produce a single RGB image of equal
horizontal and vertical resolution as the input. (B) A typical workflow
example that used maximum-intensity projection (MIP) of the input Z-stack for
a compressed representation that was passed into a 2D instance segmentation
network that predicts object masks. (C) Our proposed workflow diagram of data
in the full model. The 3D data is passed to PEN, which passes its 2D RGB
output to the 2D instance segmentation network of choice that produces the 2D
predicted elements, such as instance masks.

Training of DNNs require a large amount of annotated data that has similar
characteristics, such as the resolution, cell size and spatial distribution,
to the data of interest. Where such data is not available, augmentations can
be used to achieve satisfactory performance. The training data utilized here
consists of MDA-MB-231 cells recorded with confocal microscopy. As shown in
Figure 2A, the gray-scale images are recorded at a low axial resolution of
$\Delta$z = 10 $\mu$m, whereas the x-y plane resolution is 0.538 $\mu$m/pixel.
As shown in Figure 2B, morphological features of cells are almost completely
lost in the axial dimension. The resolution discrepancy associated with the
imaging setup, which is often the preferred choice given the photon budget,
makes it particularly desirable to perform segmentation based on information
in the x-y plane.

In order to allow segmentation even in datasets of high cell densities, we
prepare annotated training images through an augmentation strategy. First, we
experimentally obtain confocal images of low cell densities. The confocal
stacks cover an axial range of 120 $\mu$m, at steps of 10 $\mu$m, where cells
rarely overlap when projected on x-y plane. This makes it easy to annotate
cells automatically using simple contrast-based segmentation and manually
correct for errors (Figure 2C top). We then augment the data by artificially
duplicating cell images along with their annotations, and apply spatial
translation and rotation to combine with the original images. This created an
annotated dataset with three times higher cell density (Figure 2C bottom).

Figure 2: 3D image data for training. (A) MDA-MB-231 GFP cells are embedded in
3D collagen matrices and imaged with confocal microscopy at a low axial
resolution of $\Delta$z = 10$\mu$m, resulting in as few as two image slices
per cell. (B) MIPs taken over each coordinate axis, where the sub-sampling in
the axial dimension results in visible uncertainty in cell boundaries and
morphologies in the X- and Y-projections. (C) Linear depth projections of cell
image z-stacks; [top] training image of size 256 x 256 pixels shows 4 cells
distributed in 3D with no augmentations applied, and [bottom] same training
image with augmentation applied to increase local cell density. See
Supplementary S2 for depth projection information. Scale-bars = 30 $\mu$m.

As mentioned in the above, PEN is a module that is placed in front of a 2D
segmentation network. We first pair PEN with a modified CellPose network.
CellPose is a 2D U-Net architecture that predicts horizontal and vertical
flows along with probability maps for cell/background and cell edges for each
2D test image. To resolve multiple cells that overlap when projected on a 2D
plane, we modify CellPose to predict $N_{out}$ output channels, where
$N_{out}$ is a tunable integer hyperparameter we set to 3 for this work (see
also Supplementary S1). For each annotated cell, we assign its label to one of
the output channels as ground-truth. The output channel assignment is
determined by k-means classification of the axial positions of the annotated
cells in the image. Therefore, the depths of cells are monotonically but
nonlinearly mapped to the output channels (see also Supplementary S2). As
shown in Figure 3C-D, CellPose cannot distinguish cells that are overlapping
using maximum-projections as inputs but can correctly identify individual
segmentations of cells when PEN is trained in conjunction, as shown in Figure
3E-F.

Figure 3: Comparison of algorithm predictions of cell masks and bounding boxes
over two example images. (A-B) Ground-truth object outlines of expert labeled
MDA-MB-231 cells are shown in random colors over a MIP image. (C-D) Instance
segmentation and bounding box predictions made by MIP+CellPose, (E-F)
PEN+CellPose, (G-H) MIP+Mask-RCNN, and (I-J) PEN+Mask-RCNN. Each predicted
object is randomly colorized. Scale-bars = 10 $\mu$m.

As a comparison, we also pair PEN with Mask-RCNN (PEN+MaskRCNN). Mask-RCNN is
a DNN that consists in part of a Res-Net feature pyramid network which feeds
to a RPN that proposes bounding-box regions to later be segmented. Since the
RPN may propose overlapping bounding boxes, it may allow for a single 2D pixel
to belong to more than one object. We did not modify the output structures of
Mask-RCNN, but we expected that the addition of PEN to pass additional 3D
spatial information that Mask-RCNN could utilize to distinguish 2D instances.
As shown in Figure 3G-J, the addition of PEN does not qualitatively improve
the segmentation ability of Mask-RCNN compared to training with MIP inputs.
This is consistent with previous reports showing Mask-RCNN often struggles in
cases of overlapping instances [28], as proposed regions in Mask-RCNN during
inference are reduced using non-maximum suppression to prevent multiple
detections of the same instance.

To evaluate the performance of different network configurations, we
systematically compare four metrics that have been introduced previously [27].
The results are shown in Table 1. Specifically, we compute the Jaccard Index,
Precision, Recall, and a Quality metric which measures the segmentation
quality (see also Methods). First, consistent with previous reports, CellPose
outperforms Mask-RCNN on the Jaccard Index and has improved segmentation
quality [27]. Comparing the addition of PEN to each network, on a low-density
cell image dataset with > 4,000 annotated cells where fewer than 0.6% of cells
displayed any axial overlap with another cell, the training scheme of CellPose
using 2D MIP inputs (MIP+CellPose) slightly outperforms PEN+CellPose on most
metrics. However, when compared to a dataset consisting of high-density cell
images where 36.8% of cells had axial overlap with another cell, PEN+CellPose
greatly outperforms MIP+CellPose on recall, which measures the ability of the
network to detect and segment cells in an image with an intersection over
union threshold of 50%. The poor performance noted in precision is a result of
a high frequency of false-positives. Specifically, PEN+CellPose is prone to
multiple detections on the same cell, as a result of activation in multiple
channels of the output probability maps. On Mask-RCNN, the addition of PEN
slightly improves most metrics over both datasets with the consistent
exception of the segmentation quality, compared to MIP inputs. However, the
performance boost of PEN in Mask-RCNN is less appreciable in comparison to its
application in CellPose, particularly in recall on the high-density dataset.

Table 1: Quantitative performance instance segmentation networks with the Projection-Enhancement Network (PEN). Performance of CellPose, a U-Net style network, and Mask-RCNN, a region-based network, were evaluated when trained on MIPs or in conjunction with PEN. Models were evaluated on a low-density cell dataset (N = 4082) with fewer than 0.6% of cells overlapping axially at an average of just 7.2% intersection over union in the MIP, and a high-density cell dataset (N = 111) where 36.8% of cells were overlapping axially at an average 12.7% intersection over union. Metrics are measured at a minimum intersection over union of 50% for true-positive detections. See Methods for details regarding metrics. | Low-Density | High-Density
---|---|---
Model | Jaccard | Precision | Recall | Quality | Jaccard | Precision | Recall | Quality
PEN+CellPose | 0.523 | 0.574 | 0.854 | 0.807 | 0.518 | 0.616 | 0.766 | 0.782
MIP+CellPose | 0.656 | 0.729 | 0.869 | 0.853 | 0.432 | 0.731 | 0.514 | 0.727
PEN+Mask-RCNN | 0.591 | 0.700 | 0.791 | 0.744 | 0.525 | 0.875 | 0.568 | 0.673
MIP+Mask-RCNN | 0.588 | 0.731 | 0.751 | 0.759 | 0.398 | 0.750 | 0.460 | 0.700

PEN maps a 3D gray scale image stack to a compressed 2D RGB representation,
and the mapping algorithm is learned by training PEN in conjunction with a
downstream network. To understand the dependence of PEN on its paired network,
we compare the output of PEN when it has been trained with CellPose and Mask-
RCNN respectively. Figure 4B shows the output of PEN when trained in
conjunction with CellPose, evaluated on a test image stack of low cell density
as input. The coloration corresponds to an approximate object depth where
cells with lower-to-higher axial positions are mapped to red-to-green-to-blue
color channels. This verifies that PEN offers a data-driven approach to color
code depth from a 3D image. Inspired by this finding, we further compare PEN
to a linear depth embedding algorithm in Supplementary S2 and show that PEN
outperforms the linear depth embedding over all metrics analyzed in this work,
as seen in Supplementary Table 1. The depth-encoded image from PEN explains
why the modules addition to CellPose boosts its performance in the high-
density dataset in comparison to training with the max-projection images, as
shown in Figure 4A. Interestingly, we find that when trained with Mask-RCNN,
PEN does not similarly encode depth, as shown in Figure 4C. Instead, we find
that PEN acts as a low-pass filter, and only maintains activations in one
color channel. The blurry low-passed image from PEN demystifies why the
quality of segmentations in Mask-RCNN are slightly worse, but makes cells in
the image easier to detect and improves precision, recall, and the Jaccard
metric.

Figure 4: Outputs of the Projection-Enhancement Network (PEN) after successful
training. (Left) A reference MIP image of MDA-MB-231 GFP cells distributed
within a 3D image stack, (Top) the output of PEN when trained in conjunction
with CellPose, and (Bottom) the output of PEN when trained in conjunction with
Mask-RCNN. Scale-bars = 50 $\mu$m.

Following successful training with augmented data (Figure 2), we test if
PEN+CellPose can handle experimental 3D images with high cell densities and
low axial resolution. To this end, we create a sample of two cancer cell
spheroids seeded next to each other in 3D collagen matrix (Figure 5 top).
After 1 day of cell invasion into the matrix, we image the sample with an
x-y-z tile scan that covers a volume of 3020 x 1492 x 120 $\mu$m3. The
resolution in the x-y plane is 0.538 $\mu$m/pixel, and the resolution in the
axial direction is 10 $\mu$m/pixel. Visually (Figure 5 bottom), the
disseminated cells are identifiable but display significant overlap in the 2D
projection. The cells within the spheroid boundary are, however, difficult to
distinguish even by an experienced researcher. We apply the trained
PEN+CellPose model to the 3D image stack. The segmented cells are randomly
colored and plotted over the original data (gray). PEN+CellPose identified
1037 cells associated with the spheroids on the left, and 667 cells associated
with the spheroids on the right. Cells disseminated from the spheroids are
well segmented. Their elongated shape and various types of protrusions, such
as fan-shaped lamellipodia and finger-shaped filopodia, are well preserved.
Not surprisingly, the model performs poorly in regions deep within the
spheroids. Therefore, we conclude that PEN enables 2D instance segmentation
networks to quantify the 3D invasion of tumor spheroids where the imaging
covers a large volume under low axial resolution.

Figure 5: PEN+CellPose instance segmentation of a 3D system of disseminated
cells from side-by-side MDA-MB-231 GFP spheroids. (Top) A MIP image of two
MDA-MB-231 GFP spheroids separated by $\approx$1 mm that were gelled in 1.5
mg/mL collagen at 37∘C and imaged immediately. (Center) The same spheroids
were imaged after 22 hours of invasion, and an overlay of the instance
segmentation is performed by PEN+CellPose and shown on top of the gray-scale
MIP. Over 1700 unique, randomly colored, detections are shown in the lower
image. (Lower Insets) Zoomed sections of each spheroid illustrate the effect
of crowding on PEN+CellPose performance. Top and Center image scale-bars = 500
$\mu$m, inset scale-bars = 100 $\mu$m.

After illustrating the application and performance of PEN for spheroid
invasion, we investigate the importance of architectural components to the
success of PEN through an ablation study, as shown in Table 2. We first
examined the contributions of the smallest ($K=1$) and largest ($K=11$)
convolutional kernel sizes. We find that removal of either kernel does not
effect the performance after retraining compared to the base PEN+CellPose
model, indicating that the successful axial encoding seen in Figure 4B results
from the intermediate kernel sizes, in agreement with the fact that most cells
in the training set typically span several slices. However, we expect that
including the range of kernel sizes allows PEN to remain robust to new
datasets with different axial resolution.

Next, we investigate replacement of the secondary convolutional block in each
branch of PEN with a max-pooling operation over the axial dimension (Branch
Max). This alteration makes the network more shallow with fewer parameters to
learn. We observed similar performance on the low-density dataset, and a
slight decrease in all metrics on the high-density dataset compared to the
base model. We speculate that the max-pooling operation makes the network
over-reliant on the initial convolution of each branch to learn to incorporate
the axial information to the output projection.

We then explore an alternative method to combine the spatial information
learned in each branch of PEN by replacing the final convolutional block with
a max-pooling operation (Collect Max). We find that the this model slightly
outperforms the base model on recall over both datasets. Here, we choose to
keep the convolution despite the comparable performance of the max-pooling
layer to maximize the expressive ability of the module, since the pooling
operation can be learned by the convolution.

Finally, we investigate the ground-truth assignment strategy used to assign
annotated cells to the $N_{out}$ output channels discussed in Supplementary
S1. In the base PEN+CellPose model, cells are assigned to $N_{out}=3$ ground
truth output channels based on their z-position, compared to the Random GT
model where cells are randomly assigned to $N_{out}=3$ channels, and the
$N_{out}$ = 1 model where cells are assigned a single output channel. We find
that random assignment results in very poor performance across all metrics of
both datasets. Additionally, by not including multiple output channels, we
increase the performance of the network on the low-density dataset as a result
of fewer false-positives, but yields a dramatically decreased performance in
recall on the high-density dataset as the network fails to detect superposed
objects. We conclude that multiple output channels are vital to the
performance of PEN, and that an assignment strategy based on cell position
allows PEN to learn and pass axial information to the downstream network.

Table 2: Ablation study of the PEN + CellPose (base) instance segmentation network. To evaluate the effects of ablation, each model was retrained from an initialized set of random weights. We evaluated removal of the $K=1$ and $K=11$ kernel sizes, thereby removing an individual branch of PEN shown in Figure 1A. The subsequent convolution in each branch was replaced with a max-pooling operation in the axial dimension in the Branch Max model. The final convolution in PEN was replaced with a max-pooling operation in Collect Max model. Finally, the ground-truth assignment strategy to the available $N_{out}$ channels of our modified CellPose algorithm was set to randomly assign cell labels to $N_{out}=3$ channels in the Random GT model, and to a single $N_{out}=1$ channel in the $N_{out}$ = 1 model. All models were evaluated on a low-density cell dataset (N = 4082) with fewer than 0.5% of cells overlapping axially at an average of just 7.2% intersection over union in the MIP, and a high-density cell dataset (N = 111) where 36.8% of cells were overlapping axially at an average 12.7% intersection over union. Metrics are measured at a minimum intersection over union of 50% for true-positive detections. See Methods for details regarding metrics. | Low-Density | High-Density
---|---|---
Model | Jaccard | Precision | Recall | Quality | Jaccard | Precision | Recall | Quality
Base | 0.523 | 0.574 | 0.854 | 0.807 | 0.518 | 0.616 | 0.766 | 0.782
\- K = 1 | 0.502 | 0.546 | 0.864 | 8103 | 0.518 | 0.610 | 0.775 | 0.782
\- K = 11 | 0.489 | 0.527 | 0.871 | 0.812 | 0.449 | 0.512 | 0.784 | 0.762
Branch Max | 0.5202 | 0.564 | 0.870 | 0.818 | 0.4785 | 0.600 | 0.703 | 0.746
Collect Max | 0.485 | 0.520 | 0.877 | 0.816 | 0.5298 | 0.619 | 0.802 | 0.771
Random GT | 0.009 | 0.023 | 0.014 | 0.686 | 0.0759 | 0.125 | 0.162 | 0.592
$N_{out}$ = 1 | 0.6578 | 0.734 | 0.863 | 0.840 | 0.480 | 0.811 | 0.541 | 0.710

##  3 Discussion

Biomedical research routinely produces 3D image stacks that cover a large
volume but have a low axial resolution as limited by practical considerations
such as photo damaging, and temporal resolution [17, 24]. To facilitate cell
segmentation in such datasets, here we introduce the Project-Enhancement
Network (PEN). PEN is a shallow, multiscale, convolutional neural network that
encodes a 3D image stack to a 2D RGB color image, which can be subsequently
passed to a 2D segmentation algorithm, as shown in Figure 1. We show that when
paired with state-of-the-art DNNs of 2D segmentation, PEN enables accurate
detection of cells densely populated in 3D image stacks of low axial
resolutions, as illustrated in the examples in Figure 3.

In the training of PEN we take a strategy that leverages data augmentation,
which avoids tedious manual labeling to generate annotated data [25]. We find
the strategy very effective and can be easily automated by first segmenting
low density cell images, then augmenting to artificial high density images, as
in the example of Figure 2. Employing this training strategy, we show that
PEN+CellPose network can simultaneously detect over one thousand breast cancer
cells disseminating from tumor spheroids, as seen in Figure 5.

We find that the performance of PEN depends on the downstream network it is
paired with. In this work, we compared the performance of PEN in conjunction
with two leading DNNs in cell-science, CellPose and Mask-RCNN [27, 2], as
computed in Table 1. Significantly, we found that Mask-RCNN did not result in
improved performance when built with PEN. A major structural difference in
region-based CNNs compared to U-Net style networks is the extraction of
regions for segmentation, here through a region proposal network (RPN) in
Mask-RCNN. To make the algorithm more efficient, the developers of the RPN in
Mask-RCNN chose 3 size-scales and 3 aspect-ratios for the k-anchor boxes
proposed within each sliding window [2]. While the network is therefore robust
against translations, random orientations and high variance in morphology make
many cell-image datasets difficult to determine best size and aspect ratio
parameters. In contrast, the efficacy of PEN is purely data-driven and does
not restrict object orientation or scale. Furthermore, the RPN has its own
loss function to minimize, whereas PEN is only subject to the learning
objectives of the instance segmentation network it is attached to and the data
that is used as training. On one hand, no additional loss function is a
feature of PEN, making it light-weight and a plug-and-play module. However, on
the other hand, no direct learning objective makes PEN susceptible to learn
inconsistent or poor feature embeddings as a result of underlying patterns in
the data. Taken all together, we suggest PEN to be paired with non-region-
based downstream networks, specifically U-Net style segmentation networks.

Our results shed light on the explainability of DNNs [30], as visualized in
Figure 4. In the PEN+CellPose configuration, we show that after training, PEN
learns to become a nonlinear depth encoder. This makes it possible for the 2D
CellPose to detect overlapping cells on a 2D plane using the depth-encoding
color information. In the PEN+MaskRCNN configration, however, PEN learns to
become a low-pass filter. We speculate that the non-maximum suppression used
in region-proposal networks to filter out multiple detections of objects with
significant intersection-over-union prevents Mask-RCNN from detecting
overlapping cells in any 2D projected image. However, the learned embedding
helps to improve the segmentation of the single detected object, as edges are
more easily distinguished in the low-pass image. Therefore, we find that after
training, PEN turns an input image into a semantic embedding that represents
the best image transformation to maximize the learning objectives of the
neural network it is attached to.

Through a systematic ablation study in Table 2, we find the performance of
PEN+CellPose critically depends on the assignment strategy of ground-truth
annotations to multiple $N_{out}$ predicted channels. The modifications of
CellPose in this work, particularly expanding the predicted maps to multiple
channels corresponding to object depth in the 3D image stack, are vital to
detect overlapping cells in 3D. Indeed, reducing $N_{out}$ from 3 to 1
seriously deteriorates the segmentation performance. It is interesting for
future studies to further explore the optimal $N_{out}$ that balance the
computational cost and segmentation power.

In conclusion, we propose PEN as a plug-and-play module that provides a data-
driven approach to compress a 3D image stack into a 2D RGB representation as
inputs for 2D instance segmentation networks. We highlight PEN’s utility in
the detection of disseminated cells from cell-dense spheroids and in settings
of significant cell-cell overlap. Our result is a deep-learning solution for
instance segmentation in a data regime often overlooked in the field. We
envision PEN to be a useful tool for a wide range of applications such as in
research of cancer and developmental biology.

##  4 Methods

###  4.1 Maintenance of MDA-MB-231 GFP Cells

GFP-labeled MBA-MB-231 human breast carcinoma cells are purchased from
GenTarget Inc. and are maintained according to the manufacturer’s
instructions. Briefly, growth media is prepared using Dulbecco’s Modified
Eagle Medium (Gibco, US) supplemented with 10% fetal bovine serum (Gibco, US),
1% penicillinstreptomyocin (Gibco, US), and 0.1 mM non-essential amino acid
(NEAA 100x, ThermoFisher, US). Generally, cells are cultured at less than 80%
confluency and seeded on culture dishes at recommended concentrations and
maintained for up to 12 passages. Cells are kept in culture flasks in a tissue
culture incubator at 37∘C and 5% CO2.

###  4.2 3D Cell Culture

Training images were acquired from experiments of GFP-labeled MDA-MB-231 cells
dispersed in 3D collagen matrices at a low cell-density. For these
experiments, collagen solutions were prepared by diluting rat-tail collagen
type I (Corning, US) with prepared growth medium, phosphate-buffered saline
(PBS, 10x), and sodium hydroxide (NaOH, 0.1M) to a concentration of 1.5 mg/mL
or 3.0 mg/mL with pH 7.4. To embed the cells in 3D collagen matrices, cells
are suspended at very low density of approximately 650 cells/$\mu$L in ice-
cold neutralized collagen solution and added to a 35 mm collagen coated glass
bottom dish with a 7 mm microwell diameter (No. 0 coverslip, MatTek, US). The
microwell containing ice-cold cell-collagen solution is covered with a
coverslip so that the dish may be inverted during gelation to ensure
dispersion of cells in 3D. The dish is then incubated on either a warming
plate set to 25∘C, or in a tissue culture incubator (37∘C, 5% CO2) for 30
minutes in order to solidify the matrix. The coverslip is removed after
gelation time and the cellularized ECM is immersed with tissue culture medium
and continuously incubated for 24 hours before imaging. Prior to imaging, 1 M
HEPES (Gibco, U.S.) is added to 10% v/v to DMEM growth media in microwell
MatTek dishes containing cellularized ECM to maintain pH during imaging. The
dish is then imaged as described in Methods.

###  4.3 3D Spheroid Culture

GFP-labeled MDA-MB-231 spheroids seen in Figure 5 are grown following methods
by Thermo Fisher Scientific [1]. Briefly, MDA-MB-231 cells are first cultured
and seeded at low density in growth medium (100 cells/$\mu$L) in a 96-well
low-attachment dish (manufacturer). Cells are collected at the bottom of the
dish by centrifuging at 290 g for 3 minutes. After overnight culture in a
tissue culture incubator (37∘C, 5% CO2), rat-tail collagen type I (Corning,
US) is added to each well to a final concentration of 6 $\mu$g/mL in order to
promote compact spheroid growth. The 96-well plate is again centrifuged at 100
g for 3 minutes, and placed back into the tissue culture incubator. Between
3-5 days later, spheroids are cultured as follows. First, an ice-cold
neutralized collagen solution is prepared as previously described to a
collagen concentration of 1.5 mg/mL and final pH of 7.4 and kept on ice.
Spheroids are then detached from the 96 well-plate by gentle expulsion of
growth media with a pipette. Once the spheroid is visibly free-floating, the
spheroid is gently pipetted out of the well and expelled into the ice-cold
collagen solution. Multiple spheroids may be added to the same collagen
solution as desired. The spheroid-collagen solution is then added to 35 mm
collagen coated glass bottom dish with a 7 mm microwell diameter (No. 0
coverslip, MatTek, US). The dish is then placed into the tissue-culture
incubator for 15 minutes to solidify the matrix, then removed and immersed in
DMEM culture medium with 1M HEPES (Gibco, U.S.) added to 10% v/v. As invasion
proceeds rapidly (within hours), the dish is immediately taken for imaging.

###  4.4 Microscopy

3D imaging is done with a Leica TCS SPE confocal microscope with a 20x oil
immersion lens (NA 0.60) equipped with a stage-top incubator (Ibidi).
Generally, experiment dishes are placed on an on-stage incubator (Ibidi) which
maintains a constant 37∘C temperature during imaging. A drop of the immersion
oil (type HF, Carquille, U.S.) is left in contact between the dish and the
objective lens to equilibrate for an additional half hour to prevent drift
while imaging. The acquired raw images are gray-scale with a resolution of
1024 x 1024 pixel2. The voxel size has been calibrated to equal 0.538 $\mu$m.
A single x-y plane is imaged every 10 $\mu$m in the z-dimension per
experiment, resulting in as few as 2 images per cell depending on the
orientation and morphology of the cell.

###  4.5 Quantitative Metrics

In order to make direct comparisons, we have used the same analysis used in
the original CellPose work to analyze the performance of CellPose and Mask-
RCNN networks with and without the addition of PEN [27]. Briefly, predicted
objects are assigned to ground-truth labels, and thus labeled as true-
positives, using a linear sum assignment to minimize the intersection over
union loss. However, we require predicted objects to have an intersection over
union with their corresponding ground-truth assignment of 0.5 to be an
eligible candidate for a true-positive (TP) label. All non-matched predictions
are labeled as false-positives (FP), while all missed ground-truth objects are
labeled as false-negatives (FN). In this work, the metric of Jaccard Index is
the same definition as average precision used in the original CellPose work,
defined as

$Jaccard\>Index=\frac{TP}{TP+FP+FN}.$ (1)

Precision is used to measure the percentage of false-positive predictions made
by the deep-learning segmentation, and defined as

$Precision=\frac{TP}{TP+FP}.$ (2)

Recall measures the ability of the network to detect and segment objects in an
image, and is defined as

$Recall=\frac{TP}{TP+FN}.$ (3)

Finally, because the definition of Precision here does not measure the ratio
of properly identified pixels to all predicted object pixels, we include a
final metric we call "Quality". This metric measures the average segmentation
quality, and is defined as the average intersection over union of true-
positive elements, or

$Quality=\overline{IoU}_{matched}.$ (4)

###  4.6 Dataset Information

For the work reviewed here, the acquired images have not undergone any
additional image processing – a testament to the effectiveness of the deep-
learning networks to detect cells. In order to process z-stacks, image planes
are stored within OME-TIFF file formats using Tifffile Python library. Ground-
truth annotations were acquired by manual thresholding of z-stacks and taking
MIPs of individual cells. Specifically, prior to thresholding, fluorescence
images are background subtracted using a rolling ball radius of 50 pixels
(26.88 $\mu$m) and then log-transformed in order to make cell edges highly
visible and so that less fluorescent-intense cells are also quantified using
ImageJ (NIH). A manual threshold is then applied for each image. After, cells
are manually segmented for each z-stack if applicable. Since consecutive
z-stacks may have cell overlap, custom Matlab scripts are then used to
determine if the same cell is in multiple z-stacks. After, we take a MIP (2D)
of each cell. We then save the cell masks as sets of vertices using scikit-
image for a compact representation of the data in a JSON format. The ground-
truth vertices are then imported and converted back to masks and cell borders,
and horizontal and vertical gradients are calculated for the CellPose
algorithm using the heating algorithm described by [27]. A subset of the
training and validation data, the curated high cell density data, and the
spheroid image of Figure 5 is shared on Figshare at [20, 19, 9]. The full
training set is available upon request.

###  4.7 Network Training Procedures

All CellPose networks with and without PEN described in this work were trained
on a single Nvidia Tesla K80 12gB GPU. Training images were 256 x 256 pixel2
in size, and a batch size set to 8 images, and networks were trained for 50
epochs with 50 iterations per epoch. The results discussed in this work were
taken from the model weights minimizing the total loss of a 100 image
validation set, which underwent the same augmentation and cropping prescribed
for the training set. The total loss function for CellPose networks was a
summation of cross-entropy loss for the cell-background probability map, a
mean-squared error loss for the predicted horizontal and vertical gradients,
and a dice-loss for the probability map of cell edges. Training used
stochastic gradient descent with momentum set to 0.9, with a learning rate of
0.02, a weight decay of 1e${}^{-}5$ for regularization, and gradient clipping
of 5 to prevent exploding gradients during.

All Mask-RCNN networks are similarly trained on a single Nvidia Tesla K80 12gB
GPU. Training images are 512 x 512 pixel2 in size, and a batch size set to 2
images, and networks were trained for 50 epochs with 50 iterations per epoch.
We use the same weighting scheme for losses and follow the same training
procedure as used in the original implementation by [2]. We use a ResNet-50
backbone, RPN anchor scales of 8, 32, 64, 128, and 256, and anchor ratios of
0.5, 1, and 2. We use an NMS threshold of 0.9 during training, and reduce the
threshold to 0.7 during inference to consider more proposed regions.

###  4.8 Code Availability

The PEN was developed to be a simple plug-and-play module, easily implemented
on top of any 2D instance segmentation network that accepts the typical 2D RGB
input image structure. We are in the processing of developing PEN into an
installable Python library through the Python Package Index. CellPose was
developed by Stringer and Pachitariu and originally written in PyTorch [27].
We have translated the open-source code to Tensorflow/Keras and have made
several modifications, such as prediction of cell-edges, the multichannel
output discussed in Supplementary S1, and a faster post-processing flow
algorithm, but we make no claim on their intellectual property. Mask-RCNN was
developed by He and Girshik and the implementation developed with
Tensorflow/Keras by Abdullah [16, 2]. We have only modified the configuration
files to import the data structures discussed in this work, and to build PEN
on top of their 2D network. Details of training can be found in Methods. All
source code developed in this work and trained models are available at
https://github.com/eddy6081/PEN.

## References

  * [1] Mda-mb-231 cell line spheroid generation and characterization for ht assays.  https://www.thermofisher.com/ca/en/home/references/protocols/cell-culture/3-d-cell-culture-protocol/mda-mb-231-cell-line-spheroid-generation.html.
  * [2] W. Abdulla.  Mask r-cnn for object detection and instance segmentation on keras and tensorflow.  https://github.com/matterport/Mask_RCNN, 2017.
  * [3] A. Araujo, W. Norris, and J. Sim.  Computing receptive fields of convolutional neural networks.  Distill, 2019.  https://distill.pub/2019/computing-receptive-fields.
  * [4] O. Chapelle, B. Scholkopf, and A. Zien.  Semi-supervised learning (chapelle, o. et al., eds.; 2006)[book reviews].  IEEE Transactions on Neural Networks, 20(3):542–542, 2009.
  * [5] Ö. Çiçek, A. Abdulkadir, S. S. Lienkamp, T. Brox, and O. Ronneberger.  3d u-net: learning dense volumetric segmentation from sparse annotation.  In International conference on medical image computing and computer-assisted intervention, pages 424–432. Springer, 2016.
  * [6] M. de Leeuw den Bouter, G. Ippolito, T. O’Reilly, R. Remis, M. van Gijzen, and A. Webb.  Deep learning-based single image super-resolution for low-field mr brain images.  Scientific Reports, 12(1):1–10, 2022.
  * [7] R. M. Donovan-Maiye, J. M. Brown, C. K. Chan, L. Ding, C. Yan, N. Gaudreault, J. A. Theriot, M. M. Maleckar, T. A. Knijnenburg, and G. R. Johnson.  A deep generative model of 3d single-cell organization.  PLoS computational biology, 18(1):e1009155, 2022.
  * [8] E. Dusch, T. Dorval, N. Vincent, M. Wachsmuth, and A. Genovesio.  Three-dimensional point spread function model for line-scanning confocal microscope with high-aperture objective.  Journal of microscopy, 228(2):132–138, 2007.
  * [9] C. Eddy.  Mda-mb-231 spheroid z-stack, Sep 2022.
  * [10] C. Z. Eddy, H. Raposo, A. Manchanda, R. Wong, F. Li, and B. Sun.  Morphodynamics facilitate cancer cells to navigate 3d extracellular matrix.  Scientific reports, 11(1):1–10, 2021.
  * [11] A. Elhayek, M. Welk, and J. Weickert.  Simultaneous interpolation and deconvolution model for the 3-d reconstruction of cell images.  In Joint Pattern Recognition Symposium, pages 316–325. Springer, 2011.
  * [12] R. Fernandez, P. Das, V. Mirabet, E. Moscardi, J. Traas, J.-L. Verdeil, G. Malandain, and C. Godin.  Imaging plant growth in 4d: robust tissue reconstruction and lineaging at cell resolution.  Nature methods, 7(7):547–553, 2010.
  * [13] K. A. Gerbin, T. Grancharova, R. M. Donovan-Maiye, M. C. Hendershott, H. G. Anderson, J. M. Brown, J. Chen, S. Q. Dinh, J. L. Gehring, G. R. Johnson, et al.  Cell states beyond transcriptomics: integrating structural organization and gene expression in hipsc-derived cardiomyocytes.  Cell Systems, 12(6):670–687, 2021.
  * [14] S. Gupta, P. Arbeláez, R. Girshick, and J. Malik.  Aligning 3d models to rgb-d images of cluttered scenes.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 4731–4740, 2015.
  * [15] S. Gupta, R. Girshick, P. Arbeláez, and J. Malik.  Learning rich features from rgb-d images for object detection and segmentation.  In European conference on computer vision, pages 345–360. Springer, 2014.
  * [16] K. He, G. Gkioxari, P. Dollár, and R. Girshick.  Mask r-cnn.  In Proceedings of the IEEE international conference on computer vision, pages 2961–2969, 2017.
  * [17] J. Jonkman, C. M. Brown, G. D. Wright, K. I. Anderson, and A. J. North.  Tutorial: guidance for quantitative confocal microscopy.  Nature protocols, 15(5):1585–1611, 2020.
  * [18] A. Kar, M. Petit, Y. Refahi, G. Cerutti, C. Godin, and J. Traas.  Benchmarking of deep learning algorithms for 3d instance segmentation of confocal image datasets.  PLoS computational biology, 18(4):e1009879, 2022.
  * [19] A. Naylor.  Curated high density dataset, Sep 2022.
  * [20] A. Naylor.  Subset of training and validation, Sep 2022.
  * [21] C. Ounkomol, S. Seshamani, M. M. Maleckar, F. Collman, and G. R. Johnson.  Label-free prediction of three-dimensional fluorescence images from transmitted-light microscopy.  Nature methods, 15(11):917–920, 2018.
  * [22] B. Roberts, A. Haupt, A. Tucker, T. Grancharova, J. Arakaki, M. A. Fuqua, A. Nelson, C. Hookway, S. A. Ludmann, I. A. Mueller, et al.  Systematic gene tagging using crispr/cas9 in human stem cells to illuminate cell organization.  Molecular biology of the cell, 28(21):2854–2874, 2017.
  * [23] D. Sarvamangala and R. V. Kulkarni.  Convolutional neural networks in medical image understanding: a survey.  Evolutionary intelligence, pages 1–22, 2021.
  * [24] H. Schneckenburger and V. Richter.  Challenges in 3d live cell imaging.  In Photonics, volume 8, page 275. MDPI, 2021.
  * [25] C. Shorten and T. M. Khoshgoftaar.  A survey on image data augmentation for deep learning.  Journal of big data, 6(1):1–48, 2019.
  * [26] N. Silberman, D. Sontag, and R. Fergus.  Instance segmentation of indoor scenes using a coverage loss.  In European conference on computer vision, pages 616–631. Springer, 2014.
  * [27] C. Stringer, T. Wang, M. Michaelos, and M. Pachitariu.  Cellpose: a generalist algorithm for cellular segmentation.  Nature methods, 18(1):100–106, 2021.
  * [28] S. Suh, Y. Park, K. Ko, S. Yang, J. Ahn, J.-K. Shin, and S. Kim.  Weighted mask r-cnn for improving adjacent boundary segmentation.  Journal of Sensors, 2021, 2021.
  * [29] C. Szegedy, W. Liu, Y. Jia, P. Sermanet, S. Reed, D. Anguelov, D. Erhan, V. Vanhoucke, and A. Rabinovich.  Going deeper with convolutions.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 1–9, 2015.
  * [30] E. Tjoa and C. Guan.  A survey on explainable artificial intelligence (xai): Toward medical xai.  IEEE transactions on neural networks and learning systems, 32(11):4793–4813, 2020.
  * [31] A. Vaidyanathan, M. F. van der Lubbe, R. T. Leijenaar, M. van Hoof, F. Zerka, B. Miraglio, S. Primakov, A. A. Postma, T. D. Bruintjes, M. A. Bilderbeek, et al.  Deep learning for the fully automated segmentation of the inner ear on mri.  Scientific reports, 11(1):1–14, 2021.
  * [32] M. P. Viana, J. Chen, T. A. Knijnenburg, R. Vasan, C. Yan, J. E. Arakaki, M. Bailey, B. Berry, A. Borensztejn, J. M. Brown, S. Carlson, J. A. Cass, B. Chaudhuri, K. R. Cordes Metzler, M. E. Coston, Z. J. Crabtree, S. Davidson, C. M. DeLizo, S. Dhaka, S. Q. Dinh, T. P. Do, J. Domingus, R. M. Donovan-Maiye, T. J. Foster, C. L. Frick, G. Fujioka, M. A. Fuqua, J. L. Gehring, K. A. Gerbin, T. Grancharova, B. W. Gregor, L. J. Harrylock, A. Haupt, M. C. Hendershott, C. Hookway, A. R. Horwitz, C. Hughes, E. J. Isaac, G. R. Johnson, B. Kim, A. N. Leonard, W. W. Leung, J. J. Lucas, S. A. Ludmann, B. M. Lyons, H. Malik, R. McGregor, G. E. Medrash, S. L. Meharry, K. Mitcham, I. A. Mueller, T. L. Murphy-Stevens, A. Nath, A. M. Nelson, L. Paleologu, T. Alexander Popiel, M. M. Riel-Mehan, B. Roberts, L. M. Schaefbauer, M. Schwarzl, J. Sherman, S. Slaton, M. Filip Sluzewski, J. E. Smith, Y. Sul, M. J. Swain-Bowden, W. Joyce Tang, D. J. Thirstrup, D. M. Toloudis, A. P. Tucker, V. Valencia, W. Wiegraebe, T. Wijeratna, R. Yang, R. J. Zaunbrecher, A. I. for Cell Science, G. T. Johnson, R. N. Gunawardane, N. Gaudreault, J. A. Theriot, and S. M. Rafelski.  Robust integrated intracellular organization of the human ips cell: where, how much, and how variable.  bioRxiv, 2021.
  * [33] A. Wang, Q. Zhang, Y. Han, S. Megason, S. Hormoz, K. R. Mosaliganti, J. C. Lam, and V. O. Li.  A novel deep learning-based 3d cell segmentation framework for future image-based disease detection.  Scientific reports, 12(1):1–15, 2022.
  * [34] H. Wang, Y. Rivenson, Y. Jin, Z. Wei, R. Gao, H. Günaydın, L. A. Bentolila, C. Kural, and A. Ozcan.  Deep learning enables cross-modality super-resolution in fluorescence microscopy.  Nature methods, 16(1):103–110, 2019.
  * [35] W. Wang, D. Poe, Y. Yang, T. Hyatt, and J. Xing.  Epithelial-to-mesenchymal transition proceeds through directional destabilization of multidimensional attractor.  Elife, 11:e74866, 2022.
  * [36] M. Weigert, U. Schmidt, T. Boothe, A. Müller, A. Dibrov, A. Jain, B. Wilhelm, D. Schmidt, C. Broaddus, S. Culley, et al.  Content-aware image restoration: pushing the limits of fluorescence microscopy.  Nature methods, 15(12):1090–1097, 2018.
  * [37] A. Wolny, L. Cerrone, A. Vijayan, R. Tofanelli, A. V. Barro, M. Louveaux, C. Wenzl, S. Strauss, D. Wilson-Sánchez, R. Lymbouridou, et al.  Accurate and versatile 3d segmentation of plant tissues at cellular resolution.  Elife, 9:e57613, 2020.
  * [38] H. Zhang, C. Fang, X. Xie, Y. Yang, W. Mei, D. Jin, and P. Fei.  High-throughput, high-resolution deep learning microscopy based on registration-free generative adversarial network.  Biomedical optics express, 10(3):1044–1063, 2019.


[a]Anian Altherr [a]Roman Gruber

# Hadronic vacuum polarization with C* boundary conditions

Lucius Bushnaq  Isabel Campos  Marco Catillo  Alessandro Cotellucci  Madeleine
Dale  Patrick Fritzsch  Javad Komijani  Jens Lücke  Marina Krstić Marinković
Sofie Martins  Agostino Patella  Nazario Tantalo  Paola Tavella

(August 2022)

###### Abstract

We present a progress report on the calculation of the connected hadronic
contribution to the muon $g-2$ with C⋆ boundary conditions. For that purpose
we use a QCD gauge ensemble with 3+1 flavors and two QCD+QED gauge ensembles
with 1+2+1 flavors of dynamical quarks generated by the RC⋆ collaboration. We
detail the calculation of the vector mass and elaborate on both statistical
and systematic errors.

##  1 Introduction

There has been a tantalizing discrepancy between the experimental and
theoretical values for the anomalous magnetic moment of the muon [1]. On the
experimental side, two experiments at Fermilab [2] and J-PARC [3] are aiming
to further decrease the uncertainty in the next few years. Hence, progress on
the theoretical side is crucial to keep up with experiments. The main
uncertainty in the theoretical calculation stems from the hadronic vacuum
polarization (HVP) contribution. To keep up with experimental precision, a
subpercent precision is required for the HVP contribution. While there exist
several sub-percent calculations of the connected light contribution (see Ref.
[1] for an overview) at the isospin-symmetric point, the systematic
uncertainty due to neglecting isospin-breaking effects becomes crucial.
Several collaborations have started to take these effects into account for the
anomalous magnetic moment [4, 5, 6, 7, 8, 9, 10].

There exist different ways of including these effects; for a detailed review
of some of these methods see Ref. [11, 12]. Among these is
$\text{QED}_{\mathrm{L}}$ that removes the zero modes of the photonic field by
hand [13] and thus creates a non-local field theory. Alternatively,
$\text{QED}_{\mathrm{m}}$ uses a finite photon mass as an infrared regulator
[14, 15], but requires an extrapolation to zero photon mass. In
$\text{QED}_{\infty}$ the QED part is analytically calculated in the continuum
[16, 17, 18]. In this proceedings, we choose $\text{QED}_{\mathrm{C}}$ in
which the zero modes of the photon field are absent due to C⋆ boundary
conditions [19, 20, 21, 22, 23, 24], thus allowing for a local field theory
without including an additional regulator. Using these boundary conditions, we
explore the connected HVP contribution.

After a brief introduction to the HVP in Section 2, we introduce a few aspects
of C⋆ boundary conditions in Section 3. In Section 4 we present our results of
vector masses and the HVP and conclude in Section 5.

##  2 Hadronic vacuum polarization

The HVP to the muon anomalous magnetic moment can be obtained using the time-
momentum representation [25]:

$\displaystyle
a_{\mu}^{\text{HVP}}=\left(\frac{\alpha}{\pi}\right)^{2}\int_{0}^{\infty}\mathrm{d}x_{0}G(x_{0})\tilde{K}(x_{0};m_{\mu}).$
(1)

Here $\alpha$ is the electromagnetic coupling, $\tilde{K}(x_{0};m_{\mu})$ is a
kernel function defined in Ref. [26, Eq. (8)] with the muon mass $m_{\mu}$,
and $G(x_{0})$ is expressed in terms of the two-point correlation function of
the electromagnetic current $j_{k}(x)$:

$G(x_{0})=-\frac{1}{3}\sum_{k=1}^{3}\int\mathrm{d}^{3}x\expectationvalue{j_{k}(x)j_{k}(0)}.$
(2)

The integral over time $x_{0}$ in Eq. (1) can be split into two parts to take
into account the finite-time extent on the lattice and to treat the noise that
dominates the signal for large times separately [25], see Section 4.2.

##  3 C⋆ boundary conditions

At the target precision of one percent it is crucial to include QED effects in
lattice simulations, as they are expected to be of the order of one percent. A
naive implementation of QED on a finite-volume lattice with periodic boundary
conditions does not allow to simulate charged particles because states with
non-zero electric charge violate Gauss’ law. C⋆ boundary conditions [19, 21,
20, 22, 23] provide a remedy, and they do not lead to a non-local field theory
as opposed to $\mathrm{QED}_{\mathrm{L}}$. The fields obey the following
constraints in spatial direction $\hat{k}=\hat{1},\hat{2},\hat{3}$

$\displaystyle\psi_{f}(x+L\hat{k})$
$\displaystyle=\psi_{f}^{c}(x):=C^{-1}\overline{\psi}_{f}^{T}(x)$ (3)
$\displaystyle\overline{\psi}_{f}(x+L\hat{k})$
$\displaystyle=\overline{\psi}^{c}(x):=-\psi_{f}^{T}(x)C$ $\displaystyle
U_{\mu}(x+L\hat{k})$ $\displaystyle=U_{\mu}^{c}(x):=U_{\mu}(x)^{\ast}$
$\displaystyle A_{\mu}(x+L\hat{k})$
$\displaystyle=A_{\mu}^{c}(x):=-A_{\mu}(x),$

where $\psi_{f}$ and $\overline{\psi}_{f}$ are the fermionic fields of flavor
$f$, $U_{\mu}(x)\in SU(3)$ are the QCD gauge links and $e^{iA_{\mu}(x)}\in
U(1)$ are the QED gauge links; $U_{\mu}^{*}$ denotes complex conjugation and
we denote charge conjugation by the superscript $c$. The charge conjugation
matrix $C$ obeys $C\gamma_{\mu}C^{-1}=-\gamma_{\mu}^{T}$ with the Euclidean
gamma matrices $\gamma_{\mu}$. We note that the photon field $A_{\mu}(x)$ is
antiperiodic and thus does not have a zero-momentum component by construction.

###  3.1 Correlation functions

It is useful to combine $\psi_{f}$ and $\overline{\psi}_{f}$ into a spinor
doublet

$\chi_{f}(x)=\begin{pmatrix}\psi_{f}(x)\\\ \psi_{f}^{c}(x)\end{pmatrix}$ (4)

and express the action in terms of $\chi_{f}$. When using the Wilson-Dirac
formulation with a Sheikholeslami-Wohlert improvement term (see Refs. [27, 28]
for details), we find that the fermionic action can be written as

$S=-\frac{1}{2}\sum_{f=1}^{N_{f}}\chi_{f}^{T}C\sigma_{1}D\chi_{f},$ (5)

where $\sigma_{1}$ is the first Pauli matrix acting on the spinor doublet in
Eq. (4) and $N_{f}$ is the number of flavors. We note that $D$ has $24V\times
24V$ complex components as it acts on a spinor doublet ($V$ denotes the
lattice volume). Wick contractions yield (see Ref. [23, Appendix B])

$\wick{\c{\chi}_{f}(x)\c{\chi}_{f^{\prime}}(y)^{T}}=-\delta_{f,f^{\prime}}D_{f^{\prime}}^{-1}(x;y)\sigma_{1}C^{-1}.$
(6)

To evaluate the two-point function in Eq. (1) we use

$j_{\mu}(x)=\frac{1}{2}\sum_{f=1}^{N_{f}}q_{f}\chi_{f}^{T}(x)\sigma_{3}\sigma_{1}C\gamma_{\mu}\chi_{f}(x),$
(7)

where $q_{f}$ the charge of flavor $f$, $\sigma_{3}$ the third Pauli matrix
acting on the spinor doublet. When expressed in terms of $\psi_{f}$ and
$\overline{\psi}_{f}$ , $j_{\mu}(x)$ is equivalent to the familiar expression
$\sum_{f}q_{f}\overline{\psi}_{f}(x)\gamma_{\mu}\psi_{f}(x)$. Inserting Eq.
(7) into the two-point correlation function, we find

$\displaystyle\expectationvalue{j_{\mu}(x)j_{\nu}(y)}=\expectationvalue{j_{\mu}(x)j_{\nu}(y)}_{\text{conn}}+\expectationvalue{j_{\mu}(x)j_{\nu}(y)}_{\text{disc}},$

where in this proceedings we only take the connected contractions into
account111We note that there are two combinations of Wick contractions that
contribute to the connected part.:

$\expectationvalue{j_{\mu}(x)j_{\nu}(y)}_{\text{conn}}=\frac{1}{2}\sum_{f=1}^{N_{f}}q_{f}^{2}\tr[\sigma_{3}D_{f}^{-1}(x;y)\gamma_{\nu}\sigma_{3}D_{f}^{-1}(y;x)\gamma_{\mu}],$
(8)

and neglect the disconnected contribution

$\expectationvalue{j_{\mu}(x)j_{\nu}(y)}_{\text{disc}}=\left(\frac{1}{2}\sum_{f=1}^{N_{f}}q_{f}\tr\left[\sigma_{3}\gamma_{\mu}D_{f}^{-1}(x;x)\right]\right)\left(\frac{1}{2}\sum_{f=1}^{N_{f}}q_{f}\tr\left[\sigma_{3}\gamma_{\nu}D_{f}^{-1}(y;y)\right]\right).$
(9)

This expression is similar to the one obtained in periodic boundary
conditions, the only difference being that $D^{-1}$ has $24V\times 24V$
components.

###  3.2 Lattice parameters

We use configurations generated by the RC⋆ collaboration [29] with 3+1 and
1+2+1 flavors of Wilson quarks with a clover term for both $SU(3)$ and $U(1)$.
The $SU(3)$ action is $\mathcal{O}(a)$-improved in a non-perturbative way for
QCD ensembles [30], which is valid for QCD+QED ensembles up to
$\mathcal{O}(\alpha)$, and the $U(1)$ action is (so far) tree-level improved.
All configurations use C⋆ boundary conditions (see Section 3) and are
summarized in Table 1. There is one QCD ensemble (A400a00b324) and two
ensembles with dynamically generated QCD+QED fields (A360a50b324 and
A380a07b324) at two values of the fine structure constant; one is close to the
physical value and one is larger. The charged pion mass ranges between
approximately $360~{}\mathrm{MeV}$ and $400~{}\mathrm{MeV}$ and we have
$m_{\pi^{\pm}}L\approx 2.9-3.5$. We use the reference value
$\sqrt{8t_{0}}=0.415\;\text{fm}$ [31] to obtain physical units. We refer to
Ref. [29] for the details of the generation.

ensemble | A360a50b324 | A380a07b324 | A400a00b324
---|---|---|---
flavors | 1 + 2 + 1 | 1 + 2 + 1 | 3 + 1
$\beta$ | 3.24 | 3.24 | 3.24
$\alpha_{R}$ | $0.040633(80)$ | $0.007081(19)$ | $0.0$
$m_{\pi^{\pm}}[\mathrm{MeV}]$ | $358.6(3.7)$ | $383.6(4.4)$ | $398.5(4.7)$
$a[\mathrm{fm}]$ | 0.05054(27) | 0.05323(28) | 0.05393(24)
number of used configurations | 181 | 200 | 200
Table 1: Summary of used ensembles and their parameters, compare Ref. [29] for
details. All lattices have size $32^{3}\times 64$. Ensembles with 1+2+1 flavor
decomposition have degenerate down and strange quark, ensembles with 3+1
flavor decomposition have degenerate up, down and strange quark. The name of
the ensemble — say A360a50b324 — contains the approximate charged pion mass
($360~{}\mathrm{MeV}$), the bare fine-structure constant ($0.050$) and the
coupling ($3.24$). We denote the renormalized fine-structure constant by
$\alpha_{R}$ and also tabulate the number of used configurations. Note that we
have not used all available configurations, which is about 2000 for each
ensemble. We use the reference value $\sqrt{8t_{0}}=0.415\;\text{fm}$ [31] to
obtain physical units.

##  4 Preliminary results

In this section, we first examine the signal-to-noise ratio of the two-point
function and then show values for the vector masses and the HVP. Our analysis
is based on the openQ*D program package [32].

###  4.1 Signal-to-noise

Using the three ensembles in Table 1 we calculate the two-point function in
time-momentum representation $G(x_{0})$. As illustrated in Figure 1, the
relative statistical error for simulating QCD (A400a00b324) is comparable to
the error of simulating QCD+QED at physical $\alpha$ (A380a07b324). On the
other hand, simulating QCD+QED at unphysically large $\alpha$ (A360a50b324)
yields larger errors than for physical $\alpha$. This might be due to a lower
pion mass ($360~{}\mathrm{MeV}$ for A360a50b324) compared to the other
ensembles.

Figure 1: Relative error comparison. On the left-hand side, we plot the local-
local two-point function of the electromagnetic current $G(x_{0})$, on the
right-hand side we plot the relative errors. The errors for the QCD ensemble
(blue solid line) and the QCD+QED ensemble with physical $\alpha$ have
comparable relative errors, whereas the relative error for the QCD+QED
ensemble with larger $\alpha$ is slightly larger. The number of used gauge
configurations and stochastic sources are $181$ and $10$, respectively in all
three cases.

###  4.2 Mass spectroscopy

As the signal-to-noise ratio deteriorates for large times $x_{0}$, we
introduce a $x_{0,\text{cut}}$ and replace the two-point function $G(x_{0})$
by a model function for $x_{0}>x_{0,\text{cut}}$. In finite volume the two-
point function $G(x_{0})$ can be written as a sum of exponentials with
positive coefficients. We choose to only consider the leading term in the
spectral decomposition for our model function

$\displaystyle G(x_{0})_{\text{model}}=Ae^{-m_{\text{eff}}x_{0}},$ (10)

where $m_{\text{eff}}$ is the effective mass of the corresponding vector meson
ground state in lattice units and $A$ the decay amplitude.

In order to determine the coefficients $m_{\text{eff}}$ and $A$, we perform a
$\chi^{2}$-fit to the two-point correlation function in a fit range
$I_{\text{fit}}$ where excited states have decayed sufficiently, but where
there is still a clear signal. Since we are interested in the ground state, we
apply Gaussian smearing to the sink and source point to increase the overlap
with the ground state. In addition, we apply stout smearing to the gauge
fields. The effective mass $m_{\text{eff}}$ is then extracted by using a
$1$-parameter logarithmic fit to the correlator on a subset of $40$
configurations. Finally, we use the effective mass as input for the
$1$-parameter fit of the decay amplitude to the two-point function (obtained
with point sources, without smearing). In the future we plan to improve the
model function by including excited states.

Table 2 shows a compilation of the extracted masses. The statistical error is
obtained using the jackknife method and the systematic error due to choosing a
fit range for the vector mass is estimated by varying the range
$I_{\text{fit}}$. The total error in the table equals the two error
contributions added in quadrature. The mass of the ground state for the charm
contribution is not determined because the model part of the correlator has a
negligible contribution to the HVP in that case (see Figure 2).

ensemble | quark | $m_{\text{eff}}$ | $\chi_{r}^{2}$
---|---|---|---
A360a50b324 | down/strange | 0.262(7) | 1.03
up | 0.267(8) | 0.93
A380a07b324 | down/strange | 0.265(6) | 1.07
up | 0.266(4) | 0.95
A400a00b324 | up/down/strange | 0.278(7) | 0.94
Table 2: Vector meson mass for the three ensembles: The effective mass
$m_{\text{eff}}$ is obtained by performing a $\chi^{2}$-fit of the two-point
function Eq. (2) to a single exponential Eq. (10). In the last column, we
display the $\chi^{2}$ value per degree of freedom $\chi_{r}^{2}$.

###  4.3 Hadronic vacuum polarization contribution

As a preliminary check we plot the integrand of Eq. (1) in Figure 2. Here we
use a combination of conserved and local currents, i.e., $G(x_{0})$ is
obtained from
$\expectationvalue{j_{k}^{\text{conserved}}(x)j_{k}^{\text{local}}(0)}$. The
use of local currents simplifies the calculations, but requires both a
multiplicative as well as an additive renormalization constant [33]. We set
the renormalization constant due to the local current to unity. The results
for the HVP contribution of muon $g-2$ are shown in Table 3.

Figure 2: Integrand of the HVP contribution for the A380a07b324 ensemble where we use a conserved-local two-point function. Points label the actual lattice data, solid lines symbolize interpolation (for $x_{0}<x_{0,\text{cut}}$) and dashed lines stand for the model function ($x_{0}\geq x_{0,\text{cut}}$) respectively. The systematic error occurring due to using only a single exponential for the model part is significant for the up contribution (blue) and down/strange contribution (green); for the charm contribution (black) it is negligible. ensemble | flavor | $a_{\mu}^{\mathrm{HVP}}\times 10^{10}$
---|---|---
A360a50b324 | up | $309(11)$
down/strange | $77(2)$
charm | $10.62(11)$
A380a07b324 | up | $331(7)$
down/strange | $83(2)$
charm | $9.78(10)$
A400a00b324 | up/down/strange | $319(8)$
charm | $9.97(9)$
Table 3: Results for the HVP contribution. The ensembles employ C⋆ boundary
conditions in spatial directions, have a pion mass of approximately 360, 380,
and 400 MeV respectively, and are simulated on a lattice with extent
$32^{3}\times 64$ (compare Table 1). We use a combination of conserved and
local current and set the renormalization constant to unity.

###  4.4 Error contributions

We display the error contributions for the up-quark contribution of the
QCD+QED ensemble with physical $\alpha$ in Table 4. An estimate of the
statistical Monte Carlo error is obtained using jackknife. The error occurring
in the determination of the vector mass (compare Table 2) affects the
precision of the HVP contribution as does the error on the lattice spacing.
The lattice spacing has been determined in Ref. [29], which is based on the
value of the the gradient flow scale $t_{0}$ determined in Ref. [31]; the
relative error of the lattice spacing turns out to be $\frac{\Delta
a}{a}=0.53\%$ (see Table 1). We note that the gradient flow scale has only
been determined for QCD ensembles, a scale setting in full QCD+QED remains to
be worked out; see Ref. [12]. The contribution to the uncertainty on the HVP
due to error propagation of the uncertainty in scale setting is expected to be
roughly $1.8\frac{\Delta a}{a}$ (see Ref. [26, Appendix B.2]), which agrees
with our findings. In addition to statistical errors, we have systematic
errors. Due to the exploratory nature of this proceedings we do not carry out
a continuum extrapolation nor a chiral extrapolation to physical meson masses.
Neither do we assess finite-size effects. There are further systematic errors
due to the model part. To determine the vector mass and amplitude in Eq. (10)
we select a fit range and a cutoff $x_{0,\text{cut}}$. We estimate the error
due to that choice by varying the cutoff and the fit range around our chosen
value. These errors are small compared to the error that arises due to
neglecting excited states. We estimate the contribution of excited states by
bounding the correlator from above and estimating the error as the difference
between the bounds [4]. In the future we plan to reduce the error
contributions individually, see Section 5.

| variation w.r.t. | relative error
---|---|---
statistical | jackknife |  $1.21$%
err. prop. of vector mass |  $1.36$%
err. prop. of scale setting |  $0.92$%
systematic | fit range |  $0.14$%
model cutoff $x_{0,\text{cut}}$ |  $0.03$%
excited states |  $1.20$%
total |  |  $2.37$%
Table 4: Error contributions of the up contribution for ensemble A380a007b324.
The row ’jackknife’ quantifies the error due to statistical fluctuations in
the Monte Carlo simulation. The second row quantifies the error propagation of
the uncertainty on the vector mass (see Section 4.2) and the third row the
error propagation of the uncertainty on the lattice spacing. The error due to
selecting a fit range (for the vector mass and the amplitude), and a model
cutoff $x_{0,\text{cut}}$ are quantified by varying these parameters, the
error due to neglecting excited states is estimated by using a bounding
method. We do not assess continuum extrapolation nor chiral extrapolation,
neither do we estimate finite-size effects.

##  5 Conclusion and outlook

We have presented the first calculation of the connected HVP contribution to
the muon $g-2$ using C⋆ boundary conditions. Three different ensembles are
used, a QCD ensemble and two QCD+QED ensembles with different values of the
fine structure constant. The noise level for the ensemble with physical
$\alpha$ is comparable to the ensemble with QCD only, whereas for larger
$\alpha$ the noise level increases. There remain a couple of open questions
that we plan to address in the future.

First, isospin-breaking effects have not been addressed in the QCD ensemble
`A400a00b324`. There are two methods to deal with isospin-breaking effects — a
perturbative or a stochastic approach. In the perturbative approach, we expand
the QCD+QED action around the isospin-symmetric point ($\alpha=0$,
$m_{u}-m_{d}=0$) in $\alpha$ and $m_{u}-m_{d}$ and evaluate the correlation
function obtained in that way [34, 35]. On the other hand, QED fields can be
added stochastically by multiplying QCD gauge links with random $U(1)$-phases
[36, 35]. We plan to investigate both approaches.

Second, finite-size effects for the HVP can be quantified using the Hansen-
Patella method [37], where different orders of finite-volume effects are
generated by powers of $e^{-m_{\pi}L}$ with $m_{\pi}$ the pion mass, and $L$
the linear lattice extent. For periodic boundary conditions the leading-order
term is $\mathcal{O}(e^{-m_{\pi}L})$ whereas for C⋆ boundary conditions this
leading order vanishes and finite-size effects are of the order
$\mathcal{O}(e^{-\sqrt{2}m_{\pi}L})$ [38]. In addition, there are power-law
finite-size effects due to QED, which are expected to be smaller compared to
periodic boundary conditions (see for example Refs. [22] for corrections to
baryon masses).

Third, as seen in Table 4, we need to implement additional variance-reduction
methods to achieve a sub-percent precision. Firstly, we plan to increase the
number of configurations and number of sources to decrease the statistical
error and to increase the precision of the vector mass determination. In
addition, low-mode averaging [39, 40] exploits the structure of the Dirac
operator and reduces the variance that stems from the low modes. Regarding
ensemble generation, multi-level Monte Carlo [41] reduces the variance in the
correlators exponentially in the distance. We plan to consider more
sophisticated model functions than single exponentials; this will further
decrease the error due to the model part.

## Acknowledgments

We acknowledge access to Piz Daint at the Swiss National Supercomputing
Centre, Switzerland under the ETHZ’s share with the project IDs s1101, eth8
and go22. AA’s and RG’s research is funded by Schweizerischer Nationalfonds
Project No. 200021_200866. SM received funding from the European Union’s
Horizon 2020 research and innovation programme under the Marie Skłodowska-
Curie grant agreement № 813942. AC’s and JL’s research is funded by the
Deutsche Forschungsgemeinschaft Project No. 417533893/ GRK-2575 “Rethinking
Quantum Field Theory”.

## References

  * [1] T. Aoyama, N. Asmussen, M. Benayoun, J. Bijnens, T. Blum, M. Bruno et al., _The anomalous magnetic moment of the muon in the Standard Model_ , _Physics Reports_ 887 (2020) 1 [2006.04822].
  * [2] Muon $g-2$ Collaboration collaboration, _Measurement of the positive muon anomalous magnetic moment to 0.46 ppm_ , _Phys. Rev. Lett._ 126 (2021) 141801.
  * [3] M. Abe, S. Bae, G. Beer, G. Bunce, H. Choi, S. Choi et al., _A new approach for measuring the muon anomalous magnetic moment and electric dipole moment_ , _Progress of Theoretical and Experimental Physics_ 2019 (2019)  [1901.03047].
  * [4] S. Borsanyi, Z. Fodor, J.N. Guenther, C. Hoelbling, S.D. Katz, L. Lellouch et al., _Leading hadronic contribution to the muon magnetic moment from lattice QCD_ , _Nature_ 593 (2021) 51 [2002.12347].
  * [5] ETM collaboration, _Electromagnetic and strong isospin-breaking corrections to the muon $g-2$ from lattice $\mathrm{QCD}+\mathrm{QED}$_, _Phys. Rev. D_ 99 (2019) 114502.
  * [6] Fermilab Lattice, HPQCD, and MILC collaboration, _Strong-Isospin-Breaking Correction to the Muon Anomalous Magnetic Moment from Lattice QCD at the Physical Point_ , _Phys. Rev. Lett._ 120 (2018) 152001.
  * [7] Fermilab Lattice, HPQCD, and MILC Collaborations collaboration, _Hadronic-vacuum-polarization contribution to the muon’s anomalous magnetic moment from four-flavor lattice QCD_ , _Phys. Rev. D_ 101 (2020) 034512.
  * [8] RBC and UKQCD collaboration, _Calculation of the Hadronic Vacuum Polarization Contribution to the Muon Anomalous Magnetic Moment_ , _Phys. Rev. Lett._ 121 (2018) 022003.
  * [9] C. Lehner and A.S. Meyer, _Consistency of hadronic vacuum polarization between lattice qcd and the r ratio_ , _Phys. Rev. D_ 101 (2020) 074515.
  * [10] A. Risch and H. Wittig, _Leading isospin breaking effects in the hvp contribution to $a_{\mu}$ and to the running of $\alpha$_, _arXiv_ (2021)  [2112.00878].
  * [11] A. Patella, _QED corrections to hadronic observables_ , _PoS_ LATTICE2016 (2016)  [1702.03857].
  * [12] N. Tantalo, “Matching lattice QC+ED to nature.” Talk at Lattice 2022, indico.hiskp.uni-bonn.de/event/40/contributions/847/ (2022-11-30).
  * [13] S. Uno and M. Hayakawa, _QED in Finite Volume and Finite Size Scaling Effect on Electromagnetic Properties of Hadrons_ , _Progress of Theoretical Physics_ 120 (2008) .
  * [14] M.G. Endres, A. Shindler, B.C. Tiburzi and A. Walker-Loud, _Massive photons: an infrared regularization scheme for lattice QCD+QED_ , _Phys. Rev. Lett._ 117 (2015)  [1507.08916].
  * [15] M. Clark, M. Della Morte, Z. Hall, B. Hörz, A. Nicholson, A. Shindler et al., _QED with massive photons for precision physics: zero modes and first result for the hadron spectrum_ , _PoS_ LATTICE2021 (2022) 281.
  * [16] X. Feng and L. Jin, _QED self-energies from lattice QCD without power-law finite-volume errors_ , _Phys. Rev. D_ 100 (2019) 094509.
  * [17] N.H. Christ, X. Feng, L.-C. Jin and C.T. Sachrajda, _Finite-volume effects in long-distance processes with massless leptonic propagators_ , _Phys. Rev. D_ 103 (2021) 014507.
  * [18] T. Blum, N. Christ, M. Hayakawa, T. Izubuchi, L. Jin, C. Jung et al., _Using infinite-volume, continuum QED and lattice QCD for the hadronic light-by-light contribution to the muon anomalous magnetic moment_ , _Phys. Rev. D_ 96 (2017) 034515.
  * [19] A.S. Kronfeld and U.J. Wiese, _SU(N) gauge theories with C-periodic boundary conditions (I). Topological structure_ , _Nuclear Physics B_ 357 (1991) 521.
  * [20] U.J. Wiese, _C- and G-periodic QCD at finite temperature_ , _Nuclear Physics B_ 375 (1992) 45.
  * [21] A.S. Kronfeld and U.J. Wiese, _SU(N) gauge theories with C-periodic boundary conditions (II). Small-volume dynamics_ , _Nuclear Physics B_ 401 (1993) 190.
  * [22] B. Lucini, A. Patella, A. Ramos and N. Tantalo, _Charged hadrons in local finite-volume QED+QCD with C ∗ boundary conditions_, _JHEP_ 02 (2016) 076 [1509.01636].
  * [23] M. Hansen, B. Lucini, A. Patella and N. Tantalo, _Gauge invariant determination of charged hadron masses_ , _JHEP_ 2018 (2018) 146.
  * [24] L. Polley, _Boundaries for $SU(3)_{c}\times U(1)_{el}$ lattice gauge theory with a chemical potential_, _Z. Phys. C - Particles and Fields_ 59 (1993) 105.
  * [25] D. Bernecker and H.B. Meyer, _Vector correlators in lattice QCD: Methods and applications_ , _Eur. Phys. J. A_ 47 (2011) 1.
  * [26] M. Della Morte, A. Francis, V. Gülpers, G. Herdoíza, G. von Hippel, H. Horch et al., _The hadronic vacuum polarization contribution to the muon $g-2$ from lattice QCD_, _JHEP_ 2017 (2017) 20 [1705.01775].
  * [27] A. Patella, “Dirac operator.” Notes for openQxD: gitlab.com/rcstar/openQxD/-/blob/master/doc/dirac.pdf (2022-11-30).
  * [28] A. Patella, “C* boundary conditions in openQ*D code.” Notes for openQxD: gitlab.com/rcstar/openQxD/-/blob/master/doc/cstar.pdf (2022-11-30).
  * [29] L. Bushnaq, I. Campos, M. Catillo, A. Cotellucci, M. Dale, P. Fritzsch et al., _First results on QCD+QED with C* boundary conditions_ , _arXiv_ (2022)  [2209.13183].
  * [30] P. Fritzsch, R. Sommer, F. Stollenwerk and U. Wolff, _Symanzik improvement with dynamical charm: a 3+1 scheme for Wilson quarks_ , _JHEP_ 25 (2018) 1 [1805.01661].
  * [31] M. Bruno, T. Korzec and S. Schaefer, _Setting the scale for the CLS $2+1$ flavor ensembles_, _Phys. Rev. D_ 95 (2017) 074504.
  * [32] I. Campos, P. Fritzsch, M. Hansen, M.K. Marinkovic, A. Patella, A. Ramos et al., _openQ*D code: a versatile tool for QCD+QED simulations: RC* collaboration_ , _Eur. Phys. J. C_ 80 (2020) 195.
  * [33] J.C. Collins, A.V. Manohar and M.B. Wise, _Renormalization of the vector current in QED_ , _Phys. Rev. D_ 73 (2006) 105019.
  * [34] G.M.D. Divitiis, P. Dimopoulos, R. Frezzotti, V. Lubicz, G. Martinelli, R. Petronzio et al., _Isospin breaking effects due to the up-down mass difference in lattice QCD_ , _JHEP_ 124 (2012) 1.
  * [35] P. Boyle, V. Gülpers, J. Harrison, A. Jüttner, C. Lehner, A. Portelli et al., _Isospin breaking corrections to meson masses and the hadronic vacuum polarization: a comparative study_ , _JHEP_ 153 (2017) 153 [1706.05293].
  * [36] A. Duncan, E. Eichten and H. Thacker, _Electromagnetic Splittings and Light Quark Masses in Lattice QCD_ , _Phys. Rev. Lett._ 76 (1996) 3894.
  * [37] M.T. Hansen and A. Patella, _Finite-volume and thermal effects in the leading-HVP contribution to muonic $g-2$_, _JHEP_ 29 (2020) 29 [2004.03935].
  * [38] S. Martins, “Finite-Size Effects of the Hadronic Vacuum Polarization Contribution to the Muon (g - 2) with C* Boundary Conditions.” Talk at Lattice 2022, indico.hiskp.uni-bonn.de/event/40/contributions/844/ (2022-11-30).
  * [39] T. DeGrand and S. Schaefer, _Improving meson two-point functions in lattice QCD_ , _Computer Physics Communications_ 159 (2004) 185.
  * [40] G.S. Bali, S. Collins and A. Schäfer, _Effective noise reduction techniques for disconnected loops in Lattice QCD_ , _Computer Physics Communications_ 181 (2010) 1570 [0910.3970].
  * [41] M. Dalla Brida, L. Giusti, T. Harris and M. Pepe, _Multi-level Monte Carlo computation of the hadronic vacuum polarization contribution to (g - 2)_ , _Physics Letters B_ 816 (2021) 136191.