(eccv) Package eccv Warning: Package ‘hyperref’ is not loaded, but highly
recommended for camera-ready version

11institutetext: The Hebrew University of Jerusalem
11email: {Ofir.Shifman<EMAIL_ADDRESS>
# Lost in Translation: Modern Neural Networks Still Struggle With Small
Realistic Image Transformations

Ofir Shifman  Yair Weiss

###### Abstract

Deep neural networks that achieve remarkable performance in image
classification have previously been shown to be easily fooled by tiny
transformations such as a one pixel translation of the input image. In order
to address this problem, two approaches have been proposed in recent years.
The first approach suggests using huge datasets together with data
augmentation in the hope that a highly varied training set will teach the
network to learn to be invariant. The second approach suggests using
architectural modifications based on sampling theory to deal explicitly with
image translations. In this paper, we show that these approaches still fall
short in robustly handling ’natural’ image translations that simulate a subtle
change in camera orientation. Our findings reveal that a mere one-pixel
translation can result in a significant change in the predicted image
representation for approximately 40% of the test images in state-of-the-art
models (e.g. open-CLIP trained on LAION-2B or DINO-v2) , while models that are
explicitly constructed to be robust to cyclic translations can still be fooled
with 1 pixel realistic (non-cyclic) translations 11% of the time. We present
Robust Inference by Crop Selection: a simple method that can be proven to
achieve any desired level of consistency, although with a modest tradeoff with
the model’s accuracy. Importantly, we demonstrate how employing this method
reduces the ability to fool state-of-the-art models with a 1 pixel translation
to less than 5% while suffering from only a 1% drop in classification
accuracy. Additionally, we show that our method can be easy adjusted to deal
with circular shifts as well. In such case we achieve 100% robustness to
integer shifts with state-of-the-art accuracy, and with no need for any
further training.

###### Keywords:

Robustness Translation Invariance Neural Networks

##  1 Introduction

(a) DINOv2’s output score for the true label ("golf ball") as a function of
the crop location.

(b) Five Examples of different crops, the colored circles connect the image to
it’s predicted probability

Figure 1: The modern DNN DINOv2 [23] still suffers from significant changes in
the output probability for the true label ("golf ball") as a function of minor
realistic translations. (a) shows DINOV2’s probability for "golf ball" for 32
crops, each measuring 224x224 pixels, extracted from the same 256x256 pixels
image with lateral translation only. (b) Demonstrates that this translation by
a few pixels translation is nearly imperceptible. We limit our experiments to
evaluation set in which the true label object remains fully observable within
all the assessed crops, as can be seen here.

(a) First row: 3 different images cropped out from the same image and differ
from each other by 1 pixel shifts. Second row: The retrieved images that are
closest in terms of DINOv2’s representation to the query image. The search is
performed over all 1300 images in the same class. Third row: The nearest
neighbor of the top images using Open Clip representation.

(b) First & Third rows (left to right): A test image, the image translated
with cyclic shift, and with a realistic shift. One might have hard time seeing
the difference since in both cases the shift is by a single pixel. Second &
Forth rows: The images from the training set whose AFC representation is the
nearest neighbor of the top images’ representations.

Figure 2: Modern Neural Networks still suffer from inconsistency to minor
translations: (a) shows an example in which both DINO-V2 and Open-CLIP fail to
be consistent on these indiscernible translations. While (b) shows that a
method that achieves 100% consistency to cyclic shifts[21] still suffers from
lack of robustness to imperceptible realistic shifts.

Neural networks have shown remarkable performance in image classification and
have even been described as achieving "super-human" performance [12]. One
reason for the success of Convolutional Neural Networks (CNNs) is the
inductive bias that arises from the use of convolutions and pooling [8, 18]
which should in principle make them robust to small perturbations of the image
and mainly to image-translations.

Despite their reliance on convolutions, in the last five years it has been
shown repeatedly that state-of-the-art deep CNNs are not invariant to tiny
translations of the input [2, 39], and can often change their output as a
result of a simple translation of the image by a single pixel. Importantly,
recent works suggest that vision transformers also exhibit vulnerability to
such transformations[11, 29, 36]. These results add to the growing literature
on the lack of robustness of deep neural networks to small changes in the
input or to test situations that are different from the training scenario [15,
16, 27, 9, 34].

In order to make modern classifiers robust to tiny perturbations, two types of
solutions have been proposed. The first approach suggests using huge datasets
possibly in combination with self-supervised pretraining and data augmentation
assuming that the use of a large and varied dataset will cause the classifier
to learn to be robust[11, 24, 33, 26]. For example, Gunasekar [11] used data
augmentation to increase robustness to translations and argued that "In ultra-
large scale datasets, accuracy/robustness might naturally come from dataset
size itself rather than model priors." This approach was also advocated in
Radford et al. [24] which argued that "large-scale task and dataset agnostic
pre-training combined with a reorientation towards zero-shot and fewshot
benchmarking on broad evaluation suites (as advocated by Yogatama et al. [37]
and Linzen [19]) promotes the development of more robust systems" [24]. This
"robustness by large-scale" approach often goes together with the use of the
popular vision-transformer (ViT) architecture [6, 4, 22], which is usually
pre-trained on a massive dataset of up to five billion images [33, 23].

A second approach uses sampling theory to modify CNNs so that they will be
invariant to translation. In particular, this approach is based on the
observation that standard architectures may violate the Shannon-Nyquist
sampling theorem [39, 21, 41]. Several recent papers [5, 21, 30, 25] have
shown how to modify CNN architecture so that the network is always robust to
any circular (i.e. cyclic) translation of the input: when the image is
translated to the right, the rightmost column reappears as the leftmost column
in the translated image. Similar ideas have also been extended to ViT
architectures and have made it possible to make them robust to circular
translation as well [29].

In this paper, we show that both of these existing approaches to providing
robustness still fall short of attaining translation invariance to simple,
realistic perturbations. When a camera is rotated by a small angle along its
optical center, the image is indeed shifted by a small translation. But this
translation is not cyclic. Rather a small part of the scene exits the field of
view in one side and another part of the scene enters the field of view in the
other side. In image-processing terms, this corresponds to taking a large
field of view image, and performing different crops that are slightly smaller.
The bottom of figure 1(b) illustrates what we call "realistic" image
translations: each image shown on the bottom is a slightly different crop of a
larger image. All the images are translations of each other, but they are not
circular translations. Such translations have also been referred to as "crop-
shifts" [21] or "standard-shift" [30, 29].

How well do the suggested approaches to providing robustness work for these
kind of realistic translations? As can be shown in Figs. 1 and 2 they still
fall short of providing translation invariance even for tiny translations.
Figure 1(a) shows that the output of DINOv2 [23] which was trained on billions
of images in a self-supervised manner, still changes drastically when the
image is shifted by a few pixels. Fig. 2(a) shows a similar effect for both
DINOv2 and OpenCLIP, and Fig. 2(b) shows a similar failure with AFC [21], a
recent method which was designed to be invariant to cyclic shifts, but also
shows better robustness to "realistic" or "crop-shift" translations, than any
other previous works.

Modern models often serve as foundational models for diverse tasks [3, 1, 38].
Given that downstream tasks rely on these models’ computed representations,
maintaining consistency of the representation becomes vital. As can be seen in
Fig. 2, we evaluate the robustness of the feature representation by using it
to retrieve nearest neighbors from a set of 1300 images from the same class.
For the AFC method, when the images are translated by a single pixel using
cyclic shift, the networks output stays the same as the retrieved neighbor
does not change. However when it is translated using a realistic one-pixel
translation (i.e. the first column of pixels disappears in one side and a new
column of pixels is introduced in the other side), then the retrieved image
changes drastically even though the images are perceptually identical.

While these results are meant to illustrate individual failures of existing
approaches, we present here quantitative experiments that show that these
failures are surprisingly common. We also present a simple method that can be
proven to achieve any desired level of consistency, although with a modest
tradeoff with the model’s accuracy. Importantly, we demonstrate how employing
this method reduces the ability to fool state-of-the-art classifiers to less
than 5% while suffering from approximately a 1% drop in classification
accuracy. Additionally, we show that our method can be easy adjusted to deal
with circular shifts as well. In such case we achieve 100% robustness to
integer shifts, as was previously reported in [5, 21, 30, 29], but with state-
of-the-art accuracy, and with no need for any training at all.

##  2 Defining and Measuring Robustness to Translations

We use two measures to quantify the robustness of a classifier to a set of
transformations. The first one is referred to as "Consistency" [2, 5, 21] and
measures the probability that a random shift in the set of allowed
translations (e.g. a 1-pixel translation) will change the output of the
classifier. The second one, "Adversarial Robustness" was also used by Michaeli
et al. [21] and asks whether an adversary which is only allowed to replace an
image with a translation in the set of allowed translations (e.g. 1 pixel
shifts) will be able to make the classifier change its output.

Mathematically, define $T(I,\Delta)$ as the translation of image $I$ by offset
$\Delta$ and $f_{\theta}(I)$ is the output of the classifier, then consistency
is defined as

$\mathbb{E}_{I,\Delta}\left[f_{\theta}(I)==f_{\theta}(T(I,\Delta))\right]$ (1)

where the expectation is taken over images and translations $\Delta$.
Adversarial robustness if defined as

$\mathbb{E}_{I}\left[\min_{\Delta}f_{\theta}(I)==f_{\theta}(T(I,\Delta))\right]$
(2)

Here the expectation is taken over images, but the translation $\Delta$ is
chosen by the adversary in order to fool the classifier rather than randomly.
By definition, adversarial robustness is lower than consistency.

We evaluated robustness and consistency using two different definitions of the
classifier’s output $f_{\theta}(\cdot)$. In the standard ImageNet setting, the
output of the classifier is one of 1000 labels and a model is robust if this
label does not change with small translations. The problem with this
definition is that it measures the robustness of the final classification and
not the robustness of the underlying representation computed by the network,
which can be used as input for various tasks when working with DNN as a
foundation model [3, 1, 7]. Note that if the network’s representation is
provably robust, then we will always know that the classification will be
robust. However, a classifier that is robust on one task, is not guaranteed to
calculate a robust representation of an arbitrary image, especially when
dealing with data that is out of distribution [14].

Therefore, in our second definition we focus on the representation computed by
the network, which is typically a vector in dimension between 512 and 2048. A
natural approach would be to check if this embedding vector does not change
when the image is translated by a single pixel (or changes by a negligible
amount). In order to define what constitutes a negligible change in the
embedding vector, we use the vectors to retrieve the most similar image in a
dataset of 1300 images (all other images in the same class). Now the
robustness measures whether the index of the retrieved image is the same for
an original image and a translated image.

As mentioned in the introduction, we use a definition of translation that is
different from the cyclic translations used in many previous papers and
captures the effect of a small natural rotation of the camera around its
optical center. We believe this definition is more realistic and relevant to
most applications. For evaluation we follow the standard protocol in training
classifiers on ImageNet: images are resized to size $256\times 256$ and then a
crop of size $224\times 224$ is taken. An image is translated one pixel to the
right by moving the center of the crop one pixel to the right. This means that
the column of pixels at the left hand side of the image disappears from the
image, while a new column of pixels appears in the right hand side. All the
other columns are shifted to the left.

The top part of Tab. 1 shows the evaluation of different definitions of
robustness for different baseline models and different shift sizes. open-CLIP
and DINOv2 are representatives of the "robustness by large scale" approach
while AFC serves as the representative of the signal processing approach that
is guaranteed to be invariant to cyclic translations. As can be seen, all
methods are still far from robust even for one pixel translations and the best
model in terms of accuracy (DINOv2) has a 1-NN consistency of less than 85%
and an adversarial robustness less than 63% to 1-pixels shifts.

Model | Accuracy |  Adv-Rob (1-NN) Shift Size: |  Adv-Rob (Class) Shift Size: |  Consistency (1-NN) Shift Size: |  Consistency (Class) Shift Size:
---|---|---|---|---|---
|  | 1 | 3 | 5 | 9 | 1 | 3 | 5 | 9 | 1 | 3 | 5 | 9 | 1 | 3 | 5 | 9
vanila vit | 76.25* | 37.21 | 14.52 | 8.93 | 6.40 | - | - | - | - | 71.02 | 61.90 | 56.95 | 52.43 | - | - | - | -
Resnet50 | 76.74 | 46.33 | 25.08 | 17.44 | 10.39 | 86.45 | 78.16 | 74.08 | 69.78 | 75.96 | 70.23 | 67.24 | 63.44 | 94.50 | 93.13 | 92.42 | 91.60
ConvNext-T | 78.29 | 59.16 | 37.53 | 29.07 | 21.20 | 89.79 | 82.37 | 78.77 | 75.18 | 83.09 | 77.54 | 74.29 | 70.86 | 96.10 | 94.64 | 93.72 | 92.79
ConvNext-L | 83.47 | 65.58 | 45.70 | 37.35 | 29.13 | 94.35 | 89.82 | 87.62 | 85.65 | 86.33 | 82.15 | 79.75 | 77.38 | 97.92 | 97.19 | 96.76 | 96.33
open-CLIP | 69.33* | 57.01 | 35.28 | 27.81 | 22.22 | 85.37 | 74.55 | 69.40 | 65.25 | 81.98 | 76.16 | 72.99 | 70.25 | 94.44 | 92.53 | 91.47 | 90.54
DINOv2 | 84.11 | 62.42 | 40.74 | 32.68 | 26.91 | 93.45 | 88.50 | 86.26 | 84.16 | 84.77 | 79.26 | 76.62 | 74.79 | 97.60 | 96.66 | 96.21 | 95.87
AFC | 74.91 | 89.77 | 76.54 | 67.35 | 55.01 | 97.28 | 93.67 | 91.11 | 87.67 | 96.46 | 93.65 | 91.37 | 87.80 | 99.05 | 98.29 | 97.64 | 96.57
Theoretical Bound | - | 91.02 | 75.37 | 62.33 | 42.31 | 91.02 | 75.37 | 62.33 | 42.31 | 95.40 | 86.83 | 79.01 | 65.37 | 95.40 | 86.83 | 79.01 | 65.37
R-MH (ConvNext-T) | 72.87 | 94.30 | 83.93 | 74.62 | 58.12 | 97.94 | 94.01 | 90.08 | 82.95 | 97.85 | 95.15 | 92.59 | 87.80 | 99.23 | 98.25 | 97.23 | 95.21
R-R (ConvNext-T) | 73.31 | 93.47 | 81.54 | 70.84 | 52.85 | 97.70 | 93.46 | 89.66 | 82.63 | 97.52 | 94.38 | 91.39 | 85.83 | 99.14 | 98.06 | 97.09 | 95.34
R-MH (ConvNext-L) | 78.93 | 94.98 | 85.58 | 76.89 | 61.04 | 98.50 | 95.68 | 92.80 | 87.36 | 98.11 | 95.71 | 93.38 | 88.82 | 99.44 | 98.75 | 98.05 | 96.54
R-R (ConvNext-L) | 78.10 | 94.48 | 84.18 | 74.63 | 58.43 | 98.28 | 94.96 | 91.76 | 86.20 | 97.91 | 95.25 | 92.70 | 87.83 | 99.36 | 98.53 | 97.71 | 96.26
R-R (open-CLIP) | 67.77* | 95.15 | 86.26 | 77.95 | 63.38 | 97.91 | 94.06 | 90.43 | 83.75 | 97.97 | 95.39 | 92.87 | 88.19 | 99.17 | 98.14 | 97.23 | 95.64
R-MH (open-CLIP) | 65.75* | 94.63 | 84.58 | 75.22 | 59.09 | 97.79 | 93.71 | 90.07 | 83.20 | 98.17 | 95.90 | 93.65 | 89.38 | 99.21 | 98.26 | 97.34 | 95.54
R-R (DINOv2) | 82.66 | 94.73 | 85.04 | 76.21 | 61.47 | 98.68 | 96.19 | 94.01 | 89.88 | 98.00 | 95.51 | 93.13 | 88.76 | 99.51 | 98.89 | 98.34 | 97.40
R-MH (DINOv2) | 82.43 | 95.14 | 86.10 | 77.78 | 63.08 | 98.78 | 96.41 | 94.06 | 89.45 | 98.16 | 95.83 | 93.53 | 89.26 | 99.54 | 98.95 | 98.39 | 97.17

Table 1: All models accuracy, adversarial robustness (2) and consistency (1)
in terms of nearest neighbor and of class change, for shift sizes of 1,3,5 and
9 pixels. Accuracy is evaluated with linear classification (*besides open-CLIP
and vanila-ViT that are evaluated with K-NN) on a random sample of 1000 images
from ImageNet validation set. R-R stands for RICS-rand and R-MH for the use of
Mexican-Hat kernel.

##  3 Our Method

###  3.1 Robust Inference by Crop Selection

Figure 3: Our proposed method: Robust Inference by Crop Selection (RICS):
First assign scores to all of the different crops of the image with some
deterministic function, then choose the crop with the maximal score (only if
it’s also a local maximum: it has 8 neighbors, and it’s score is the highest
among them). This crop is passed through a feature extractor. Left: An image
and it’s translated version, with few crops marked. Center: The crops’ scores
as a heatmap. local maximums are highlighted with a star. Right: The chosen
crop is passed to the feature extractor.

We present a method that allows us to convert any classifier into one that is
robust to realistic image translations without the need for retraining. The
method is based on the classical object recognition paradigm that is standard
in computer vision: given an image it is customary to first find the bounding
box of the object and then pass the subimage within the bounding box to a
classifier (e.g. [10]). The problem of image classification is different from
that of object detection (e.g. many ImageNet categories are not based on a
single object on a uniform background and there is no well-defined object to
detect) but we use it as inspiration for our method.

During training of image classifiers, it is common practice to perform data
augmentation by using different crops of the same image. Our method is built
on the simple observation that we can make a classifier robust to translations
by choosing a subcrop of the image during inference, not training.

Our method, which we call Robust Inference by Crop Selection is illustrated in
Fig. 3 and Algorithm 1: Given an image, we first assign scores to all of the
different crops of the image and then choose the crop with the maximal score
(among the local maxima, i.e. a score must be larger than the score of its 8
neighboring crops and not at the edge). This crop is then passed through the
standard neural network to serve as the feature extractor of the image, and
the output of our method on the image is defined as the output to the specific
crop.

At first glance, our method seems to require a complicated score function that
will include exactly the relevant parts of the image. But in fact, in order to
ensure robustness, all we require is that the crop be chosen in a consistent
manner. Consider again figure 3 and two images of the bee that differ by a
single pixel translation. In order for our method to yield consistent
classification, we only require that the same crop be chosen for the two
images but we do not require that the crop be exactly centered on the bee.

In our implementation we use two score functions. One which we call RICS-rand
(R-R) involves performing a dot product of the crop with a fixed, random,
filter (each pixel in the filter is chosen IID from a Gaussian, and the filter
is constant for all images and crops). The output of the dot-product is then
passed through a modulo function. We use this score function to illustrate the
fact that our score function is not required to be a good object detector. A
dot product with a random filter is certainly not a good method to choose the
crop most likely to contain the object at the center, but it is still
deterministic: the score only depends on the pixel values within a crop and as
we show in our theoretical analysis, a deterministic score function that is
pseudo-random is enough to guarantee robustness of our method.

The second function is more closely related to classical methods for object
recognition but still very simplistic. We perform a dot-product with a
"Mexican Hat" filter [28, 32] and after calculating the scores of all patches,
return the crop that has the maximal score, and is also a local maximum (has a
larger score than that of the 8 neighbor crops). We start with a Mexican Hat
kernel with high standard deviation (50 pixels), that should detect big region
of interests, and if this returns scores that are very smooth, with no local
maximum (besides the edge ones), we use a smaller kernel (20 or even 9
pixels). We call this method RICS-MH (R-MH). A formal definition of our method
follows.

Algorithm 1 Robust Inference by Crop Selection

1:Input: An image $\bm{I\in\mathbb{R}^{n\times n}}$, crop size $\bm{k\times
k}$, scoring function $\bm{g:\mathbb{R}^{k\times k}\rightarrow\mathbb{R}}$,
and a classifier $\bm{f_{\theta}:\mathbb{R}^{k\times
k}\rightarrow\mathbb{R}^{d}}$

2:Output: Feature representation, $\bm{y\in\mathbb{R}^{d}}$.

3:let $G\in R^{{(n-k+1)}^{2}}$ be a matrix of the crop scores
$G_{i,j}=g(C_{i,j})$ for every crop $0\leq i,j\leq n-k$.

4:perform Non Maximum Suppression on $G$

5:let $C=\arg\max_{i,j}G_{i,j}$ be the crop who’s score is the maximal.

6:return $f_{\theta}(C)$

###  3.2 Theory

As mentioned, our method is based on the observation that we do not need to
find the "best" crop within an image before running the network. We now show
that even a random score function will guarantee a high degree of robustness.

###### Definition 1

(Pseudo-Random Crop Function): a deterministic score function for a crop is
called a "Psuedo-random crop function" if for any image, the maximum of the
score function within this image is equally likely to be at any crop of the
image.

If the crop score function were simply chosen randomly for each crop, then it
would trivially satisfy the requirement that the maximum is equally likely to
be at any crop. But we also require the function to be a deterministic
function of the pixels in the crop (similar to random hash functions).
Assuming we can find such a crop function, we can guarantee robustness.

Figure 4: Illustration of the proof of Theorem 1. If the image is $n\times n$
pixels, and the crop is $k\times k$, than from a total of $(n-k+3)^{2}$
different crops that correspond to 1 pixel translations, $(n-k-1)^{2}$ are
shared, meaning the probability that the algorithm predicts the same crop for
the image and all of its translations is
$\left(\frac{n-k-1}{n-k+3}\right)^{2}$. On the left side you can see the image
$I$ (in blue) and two versions of it translated by a single pixel colored in
yellow and red. $C_{i,j}$ is the crop that gets the highest score and returned
by algorithm 1.

###### Theorem 3.1

Consider a base classifier that classifies $n\times n$ images. Suppose that we
use the RICS method (algorithm 1) with $k\times k$ crops and a pseudo-random
score function to convert the base classifier into a new classifier, then the
probability that an adversary can alter the output of the new classifier by a
realistic translation of size at most $\Delta$ pixels is bounded by:

$p(\text{"adversary\\_succeeds"})\leq
1-\left(\frac{n-k+1-2\Delta}{n-k+1+2\Delta}\right)^{2}$ (3)

###### Proof

The proof is based on a counting argument and is illustrated in figure 4. If
the crop that is selected in algorithm 1 is the same for all translations of
the image, then the adversary will surely fail. Since algorithm 1 chooses a
crop by choosing the maximum of a pseudo-random score function, the maximum is
distributed uniformly over the image. In order for all translations to choose
the same crop, the maximum should not be at the edge of the image, as shown in
figure 4 and the probability that this happens is given by the ratio of the
areas of the center and the area of the full image.

For concreteness, suppose we are classifying $256\times 256$ images with a
crop size of $150\times 150$, then theorem 1 guarantees that any classifier of
$150\times 150$ images can be converted into a classifier whose adversarial
robustness to 1 pixel translations is at least $93\%$ using the RICS method.
Note that the result is for realistic translations: we are not assuming
cyclical translations as in [5, 21] but rather our theoretical result holds
for translations where information exits and enters the field of view. Note
also that even though our theorem is stated in terms of adversarial accuracy,
a similar result can be proven for consistency as well: any classifier can be
converted into a method that satisfies a rigorous lower bound on the
consistency with realistic transformations, and from similar counting argument
we deduce a lower bound on the consistency:

$p("consistent\\_prediction")\geq\left(\frac{n-k+1-\Delta}{n-k+1+\Delta}\right)^{2}$
(4)

Figure 5(a) shows the lower bound on the adversarial robustness to
translations of any classifier when using our RICS method as a function of the
size of the crop $k$, starting with an image of size $256\times 256$ pixels.
Note that this is a lower bound, and depending on the base classifier, the
actual robustness may be much higher. This is because even if the adversary
succeeded in changing the crop that is given to the classifier by algorithm 1,
a good classifier will often output the same output to different crops and
thus the adversary will fail to modify the classifier’s output.

(a) RICS lower bound for the adversarial robustness

(b) Accuracy of DINOv2 as a function of the crop size

Figure 5: (a) RICS lower bound for the adversarial robustness when using
images of size $224\times 224$ and shift size of $\Delta=1$. As the crop size
increases, robustness goes down. (b)The accuracy of DINOv2 for different crop
sizes. Accuracy goes up as a the crop size increases and more relevant
information exist in the image, and center crops are only slightly better than
random crops.

Figure 5(b) shows that accuracy increases as the crop size increases and that
random crops are only slightly worse than center crops. The fact that accuracy
increases but robustness decreases with crop size means that there is a
tradeoff between accuracy and robustness. In practice, different applications
will choose different crop sizes depending on how they want to address this
tradeoff.

###  3.3 Cyclic Translations

As mentioned in Sec. 1, previous works have already achieved 100% robustness
to cyclic translations using polyphase components [5, 30, 29] or polynomial
activations [21]. Despite our focus on realistic (non-cyclic) translations,
our method can be easily adapted to cyclic shifts as well, by including cyclic
crops, ensuring 100% robustness as well. The fact that our method require no
further training, and can be applied as a pre-process to a state-of-the-art
model, allows us to achieve state-of-the-art accuracy as well as 100%
robustness to cyclic-translation. The only change in the algorithm is that we
replace the definition of crops in Algorithm 1 with cyclic crops:

###### Definition 2

A cyclic crop of size $k\times k$ at location $i,j$ of an image $I$ is defined
as the top left subimage of $I$ after it has been cyclically translated so
that location $i,j$ is at the top left.

###### Theorem 3.2

Consider a base classifier that classifies $n\times n$ images. For any crop
size $k$ and any deterministic scoring function, consider using algorithm 1
with cyclic crops to convert the base classifier into a new classifier. Then
an adversary can not alter the output of the new classifier by any integer
cyclic translation. In other words, the adversarial robustness (and therefore
also the consistency) of the new classifier is 100%.

###### Proof

For any crop size $k$, the set of cyclic crops of an image $I$ and the set of
cyclic crops of the translated version $T_{cyclic}(I,\Delta)$ is exactly the
same regardless of the translation size $\Delta$ (as long as the translation
is in full pixel values, and not sub-pixel translation), and consist of
$n^{2}$ crops. Since we evaluate the score of all possible crops and choose
deterministically, then Algorithm 1 will always choose the same cyclic crop,
and hence the prediction of the classifier will be the same.

Note that introducing cyclic crops means that a crop chosen by the algorithm
may be one that is unnatural: it may include for example pixels from the
bottom of the image that become adjacent to pixels from the top of the image
and in principle, our method may decrease the accuracy of the base classifier.
Nevertheless, for any classifier, our method is guaranteed to return a new
classifier that is perfectly robust to cyclic translations. In the next
section we evaluate how much accuracy is decreased as a result of using cyclic
crops.

Figure 6: The accuracy-robustness trade-off. The left plot shows that our
method achieves the best robustness compared to exiting methods, while barely
losing accuracy comparing to the base model we use, DINOv2. The right plot
illustrates that our method can be used to improve the robustness of any model
without retraining. The adversarial robustness is based on change of the first
Nearest neighbor with 1 pixel shifts.

##  4 Results

(a) Adversarial robustness as the change of the 1-NN

(b) Adversarial robustness as the change of the class

(c) Consistency of the nearest neighbor

(d) Consistency of the predicted class

Figure 7: Adversarial robustness 2 and consistency 1 as functions of the shift
size. As we allow bigger shifts, all models tend to fail more. Adversarial
robustness is depicted as a function of the shift size. Subplot (a)
illustrates that not only our method is more robust than other methods, even
the lower bound is. Subplot (b) shows that when measuring robustness in terms
of label change, model accuracy contributes to its robustness. However, the
overall ranking of the models remains mostly unchanged. Subplots (c) and (d)
shows similar trends on the consistency metric.

We compare our method to two of the best models that have been published as
open code with trained weights: Open-Clip [33] (trained on 2 billion image-
text pairs) and DINOv2 [23], as well as AFC [21], the latest version of a
model that is constructed to be translation invariant to cyclic shifts. We
also compare to three baseline methods: ConvNext [20], ResNet50 [13] and
vanila-ViT [35] trained on ImageNet in the regular supervised manner.

For the evaluation with realistic translations we use a subset of ImageNet
consisting of 1000 images, chosen randomly from all images whose object is
known to be in all of the $224\times 224$ pixel crops (after resizing the
image to $256\times 256$). We used these images because we wanted to ensure
that all the images that we create by small translations still include the
main object and to avoid cases where the some translations push the object out
of the image. [31]

In order to enable comparison to previous work, we apply our method as if the
input is a $224\times 224$ image. For the realistic translations, we choose
the crop size $k$ in our method to be $140\times 140$ because we have found
that this leads to a very small decrease in accuracy while giving high
robustness (as can bee seen in Figure 5). In practice, ImageNet images are
usually resized to $256\times 256$ and then a $224\times 224$ crop is
selected. If our method was used to select $224\times 224$ crops, then we
would expect no decrease in accuracy compared to standard classifiers.

Fig. 6 shows that our method can be applied to any classifier, but our best
results are obtained with a baseline DINO-v2. We use both a random filter and
a Mexican hat filter. As can be seen in Fig. 6 and Tab. 1, our method is the
only one to achieve top-1 accuracy of over 80% with an adversarial robustness
of about 95% to realistic 1 pixel translations. It can also be seen that for
any base classifier, our method significantly improves the robustness with a
negligible decrease in accuracy. Figs. 7(a) and 7(b) shows the adversarial
robustness of different models as a function of the maximal allowed
translation along with the lower bound predicted by our analysis, and Figs.
7(d) and 7(c) shows the same in terms of consistency. For all range of
translations, our method is the most robust. This is true both when
adversarial robustness and consistency are evaluated using the index of the
first nearest neighbor or when they are evaluated using the class.

Cyclic Results: As stated in Theorem 2, our adjusted method achieves 100%
robustness to cyclic translations. Tab. 2 shows that applying our method our
method as preprocess barely influences Dinov2 accuracy, which is 83.9% using
RICS-MH and 84.1% without it in our experiments. These results shows that our
method has higher accuracy than all previous methods solving the 100%
consistency task.

Method | Base Model | Accuracy (%) | Consistency (%)
---|---|---|---
APS [5] | ResNet-18 | 67.6 | 100
LPS [30] | ResNet-18 | 69.11 | 100
A-SwinV2-T [29] | MViTv2 | 79.91 | 99.98
A-MViTv2-T [29] | MViTv2 | 77.46 | 100
AFC [21] | Convnext-tiny | 82.11 | 100
RICS-MH | Convnext-Large | 81.83 | 100
RICS-rand | DINO-v2 | 83.77 | 100
RICS-MH | DINO-v2 | 83.93 | 100

Table 2: Accuracy and Consistency to integer shifts for cyclic shits of our
method verses four previus works. All models are evaluated on the ImageNet
validation set (50K images), The results of the previous works in this table
are copied from the papers. All methods achieves 100% consistency, but our is
the only one that can use any pretrained model with requires no further
training, and using DINOv2 achieves the best accuracy while maintaining 100%
consistency.

Additional results and figures are attached in the appendices.

##  5 Limitations and Extensions

The main limitation of our theoretical analysis and our proposed algorithm is
that we can only guarantee robustness to small realistic image translations:
the bound decreases with the size of the allowed translation. In contrast,
when dealing with cyclical translations our method, as well as previous
methods, can guarantee 100% robustness to any circular translation. In order
to extend our theory and our algorithm to larger translations, we need to
tackle the difficult problem of going beyond a pseudo-random score function
and towards score functions that are likely to capture the meaningful parts of
the image.

One possible approach towards meaningful score functions is to note that the
RICS method with a score function that is based on dot-products can be
differentiated with respect to the kernel. We can therefore try to improve our
method by introducing a consistency loss that is applied on pairs of images
and its translations. Note this this method would still not be computationally
demanding and require only a single dot-product per crop.

Another limitation of our method is that it can only handle integer
translations. We believe that it can also be extended to fractional
translations by appropriately upsampling the image prior to computing the
crops, but we leave this investigation for future work.

##  6 Discussion

The performance of deep neural networks in image classification has been
continuously improving and one major factor allowing for this improvement has
been the use of larger and larger datasets. This has led to the belief that
the major failure modes of deep classifiers (i.e. their brittleness and
inability to generalize to scenarios that are very different from the training
set) will mostly be solved by using even larger datasets or more
augmentations. In this paper, we have focused on one example of brittleness:
the inability of modern classifier to maintain consistency when subjected to a
realistic translation of a single pixel. We have shown that methods that use
sampling theory to provably provide robustness to cyclic shifts, don’t solve
the problem of realistic shifts (where a pixel exits the field of view,
another enters on the other side). We have also shown that even classifiers
that are trained on billions of training images do not learn to be robust to
these trivial transformations. Instead of relying on sampling theory or on
large-scale datasets, we suggested a simple method that to convert any
classifier into a classifier that achieves a high degree of robustness. We
presented a theory that proves a lower bound on the robustness and this bound
by itself can be used to convert any classifier into one that has at least 93%
robustness to single pixel shifts (and even higher, with a price of accuracy),
and 100% robustness to any cyclic shift. We also presented experiments that
show that the actual performance is even better. When using DINO2 as our
baseline classifier we presented the first classifier that achieves over 82%
accuracy on ImageNet and over 95% robustness to realistic single pixel shifts
and a classifier that achieves about 84% accuracy on ImageNet and over 100%
robustness to cyclic shifts.

## References

  * [1] Awais, M., Naseer, M., Khan, S., Anwer, R.M., Cholakkal, H., Shah, M., Yang, M.H., Khan, F.S.: Foundational models defining a new era in vision: A survey and outlook. arXiv preprint arXiv:2307.13721 (2023)
  * [2] Azulay, A., Weiss, Y.: Why do deep convolutional networks generalize so poorly to small image transformations? arXiv preprint arXiv:1805.12177 (2018)
  * [3] Bommasani, R., Hudson, D.A., Adeli, E., Altman, R., Arora, S., von Arx, S., Bernstein, M.S., Bohg, J., Bosselut, A., Brunskill, E., et al.: On the opportunities and risks of foundation models. arXiv preprint arXiv:2108.07258 (2021)
  * [4] Caron, M., Touvron, H., Misra, I., Jégou, H., Mairal, J., Bojanowski, P., Joulin, A.: Emerging properties in self-supervised vision transformers. In: Proceedings of the IEEE/CVF international conference on computer vision. pp. 9650–9660 (2021)
  * [5] Chaman, A., Dokmanic, I.: Truly shift-invariant convolutional neural networks. In: Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition. pp. 3773–3783 (2021)
  * [6] Dosovitskiy, A., Beyer, L., Kolesnikov, A., Weissenborn, D., Zhai, X., Unterthiner, T., Dehghani, M., Minderer, M., Heigold, G., Gelly, S., et al.: An image is worth 16x16 words: Transformers for image recognition at scale. arXiv preprint arXiv:2010.11929 (2020)
  * [7] Fang, Y., Wang, W., Xie, B., Sun, Q., Wu, L., Wang, X., Huang, T., Wang, X., Cao, Y.: Eva: Exploring the limits of masked visual representation learning at scale. In: Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition. pp. 19358–19369 (2023)
  * [8] Fukushima, K., Miyake, S., Ito, T.: Neocognitron: A neural network model for a mechanism of visual pattern recognition. IEEE Transactions on Systems, Man, and Cybernetics SMC-13, 826–834 (1983)
  * [9] Geirhos, R., Jacobsen, J.H., Michaelis, C., Zemel, R., Brendel, W., Bethge, M., Wichmann, F.A.: Shortcut learning in deep neural networks. Nature Machine Intelligence 2(11), 665–673 (2020)
  * [10] Girshick, R., Donahue, J., Darrell, T., Malik, J.: Rich feature hierarchies for accurate object detection and semantic segmentation (2014)
  * [11] Gunasekar, S.: Generalization to translation shifts: a study in architectures and augmentations. arXiv preprint arXiv:2207.02349 (2022)
  * [12] He, K., Zhang, X., Ren, S., Sun, J.: Delving deep into rectifiers: Surpassing human-level performance on imagenet classification. In: Proceedings of the IEEE international conference on computer vision. pp. 1026–1034 (2015)
  * [13] He, K., Zhang, X., Ren, S., Sun, J.: Deep residual learning for image recognition. In: Proceedings of the IEEE conference on computer vision and pattern recognition. pp. 770–778 (2016)
  * [14] Hendrycks, D., Basart, S., Mu, N., Kadavath, S., Wang, F., Dorundo, E., Desai, R., Zhu, T., Parajuli, S., Guo, M., et al.: The many faces of robustness: A critical analysis of out-of-distribution generalization. In: Proceedings of the IEEE/CVF International Conference on Computer Vision. pp. 8340–8349 (2021)
  * [15] Hendrycks, D., Dietterich, T.: Benchmarking neural network robustness to common corruptions and perturbations. arXiv preprint arXiv:1903.12261 (2019)
  * [16] Hendrycks, D., Gimpel, K.: A baseline for detecting misclassified and out-of-distribution examples in neural networks. arXiv preprint arXiv:1610.02136 (2016)
  * [17] Krizhevsky, A., Sutskever, I., Hinton, G.E.: Imagenet classification with deep convolutional neural networks. Advances in neural information processing systems 25 (2012)
  * [18] LeCun, Y., Boser, B., Denker, J.S., Henderson, D., Howard, R.E., Hubbard, W., Jackel, L.D.: Backpropagation applied to handwritten zip code recognition. Neural computation 1(4), 541–551 (1989)
  * [19] Linzen, T.: How can we accelerate progress towards human-like linguistic generalization? arXiv preprint arXiv:2005.00955 (2020)
  * [20] Liu, Z., Mao, H., Wu, C.Y., Feichtenhofer, C., Darrell, T., Xie, S.: A convnet for the 2020s. In: Proceedings of the IEEE/CVF conference on computer vision and pattern recognition. pp. 11976–11986 (2022)
  * [21] Michaeli, H., Michaeli, T., Soudry, D.: Alias-free convnets: Fractional shift invariance via polynomial activations (2023)
  * [22] Naseer, M.M., Ranasinghe, K., Khan, S.H., Hayat, M., Shahbaz Khan, F., Yang, M.H.: Intriguing properties of vision transformers. Advances in Neural Information Processing Systems 34, 23296–23308 (2021)
  * [23] Oquab, M., Darcet, T., Moutakanni, T., Vo, H., Szafraniec, M., Khalidov, V., Fernandez, P., Haziza, D., Massa, F., El-Nouby, A., et al.: Dinov2: Learning robust visual features without supervision. arXiv preprint arXiv:2304.07193 (2023)
  * [24] Radford, A., Kim, J.W., Hallacy, C., Ramesh, A., Goh, G., Agarwal, S., Sastry, G., Askell, A., Mishkin, P., Clark, J., et al.: Learning transferable visual models from natural language supervision. In: International conference on machine learning. pp. 8748–8763. PMLR (2021)
  * [25] Rahman, M.A., Yeh, R.A.: Truly scale-equivariant deep nets with fourier layers. arXiv preprint arXiv:2311.02922 (2023)
  * [26] Rebuffi, S.A., Gowal, S., Calian, D.A., Stimberg, F., Wiles, O., Mann, T.A.: Data augmentation can improve robustness. Advances in Neural Information Processing Systems 34, 29935–29948 (2021)
  * [27] Recht, B., Roelofs, R., Schmidt, L., Shankar, V.: Do imagenet classifiers generalize to imagenet? In: International conference on machine learning. pp. 5389–5400. PMLR (2019)
  * [28] Reddy, V.K., Siramoju, K.K., Sircar, P.: Object detection by 2-d continuous wavelet transform. In: 2014 International Conference on Computational Science and Computational Intelligence. vol. 1, pp. 162–167. IEEE (2014)
  * [29] Rojas-Gomez, R.A., Lim, T.Y., Do, M.N., Yeh, R.A.: Making vision transformers truly shift-equivariant. arXiv preprint arXiv:2305.16316 (2023)
  * [30] Rojas-Gomez, R.A., Lim, T.Y., Schwing, A., Do, M., Yeh, R.A.: Learnable polyphase sampling for shift invariant and equivariant convolutional networks. Advances in Neural Information Processing Systems 35, 35755–35768 (2022)
  * [31] Russakovsky, O., Deng, J., Su, H., Krause, J., Satheesh, S., Ma, S., Huang, Z., Karpathy, A., Khosla, A., Bernstein, M., Berg, A.C., Fei-Fei, L.: ImageNet Large Scale Visual Recognition Challenge. International Journal of Computer Vision (IJCV) 115(3), 211–252 (2015). https://doi.org/10.1007/s11263-015-0816-y
  * [32] Sarvaiya, J.N., Patnaik, S.: Automatic image registration using mexican hat wavelet, invariant moment, and radon transform. IJACSA) International Journal of Advanced Computer Science and Applications, Special Issue on Image Processing and Analysis (2011)
  * [33] Schuhmann, C., Beaumont, R., Vencu, R., Gordon, C., Wightman, R., Cherti, M., Coombes, T., Katta, A., Mullis, C., Wortsman, M., et al.: Laion-5b: An open large-scale dataset for training next generation image-text models. Advances in Neural Information Processing Systems 35, 25278–25294 (2022)
  * [34] Wortsman, M., Ilharco, G., Kim, J.W., Li, M., Kornblith, S., Roelofs, R., Lopes, R.G., Hajishirzi, H., Farhadi, A., Namkoong, H., et al.: Robust fine-tuning of zero-shot models. In: Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition. pp. 7959–7971 (2022)
  * [35] Wu, B., Xu, C., Dai, X., Wan, A., Zhang, P., Yan, Z., Tomizuka, M., Gonzalez, J., Keutzer, K., Vajda, P.: Visual transformers: Token-based image representation and processing for computer vision (2020)
  * [36] Wu, H., Xiao, B., Codella, N., Liu, M., Dai, X., Yuan, L., Zhang, L.: Cvt: Introducing convolutions to vision transformers. In: Proceedings of the IEEE/CVF international conference on computer vision. pp. 22–31 (2021)
  * [37] Yogatama, D., d’Autume, C.d.M., Connor, J., Kocisky, T., Chrzanowski, M., Kong, L., Lazaridou, A., Ling, W., Yu, L., Dyer, C., et al.: Learning and evaluating general linguistic intelligence. arXiv preprint arXiv:1901.11373 (2019)
  * [38] Yuan, L., Chen, D., Chen, Y.L., Codella, N., Dai, X., Gao, J., Hu, H., Huang, X., Li, B., Li, C., et al.: Florence: A new foundation model for computer vision. arXiv preprint arXiv:2111.11432 (2021)
  * [39] Zhang, R.: Making convolutional networks shift-invariant again (2019)
  * [40] Zhou, B., Lapedriza, A., Khosla, A., Oliva, A., Torralba, A.: Places: A 10 million image database for scene recognition. IEEE transactions on pattern analysis and machine intelligence 40(6), 1452–1464 (2017)
  * [41] Zou, X., Xiao, F., Yu, Z., Li, Y., Lee, Y.J.: Delving deeper into anti-aliasing in convnets. International Journal of Computer Vision 131(1), 67–81 (2023)

## Lost in Translation - Appendixes

### Appendix 1 - More examples of DINOv2 change in predicted label

(a) DINOv2’s output score for the true label ("black widow") as a function of
the crop location.

(b) Five examples of different crops, the colored circles connect the image to
it’s predicted probability in figure (a)

Figure 8: more examples of the paper’s figure 1 showing that DINOv2 suffers
from big changes in the output probability for the true label ("black widow")
as a function of minor realistic translations. One can also note the clear
osculations in the output probability, which occurs as a result of the use of
ViT, with non overlapping patches. The translation here is horizontal.

(a) DINOv2’s output score for the true label ("water bottle") as a function of
the crop location.

(b) Five examples of different crops, the colored circles connect the image to
it’s predicted probability in figure (a)

Figure 9: more examples of the paper’s figure 1 showing that DINOv2 suffers
from big changes in the output probability for the true label ("water bottle")
as a function of minor realistic translations. One can also note the clear
osculations in the output probability, which occurs as a result of the use of
ViT, with non overlapping patches. The translation here is horizontal.

(a) DINOv2’s output score for the true label ("Yorkshire terrier") as a
function of the crop location.

(b) Five examples of different crops, the colored circles connect the image to
it’s predicted probability in figure (a)

Figure 10: more examples of the paper’s figure 1 showing that DINOv2 suffers
from big changes in the output probability for the true label ("Yorkshire
terrier") as a function of minor realistic translations. One can also note the
clear osculations in the output probability, which occurs as a result of the
use of ViT, with non overlapping patches. The translation here is vertical.

### Appendix 2 - More examples and analyses of DINOv2 change in predicted
Nearest neighbor

(a) Sixteen out of the possible $32\times 32=1024$ images of size $224\times
224$ pixels cropped from a $256\times 256$ image. The main object (lens cap)
is present in all the images.

(b) The index of the image in the training set, which is the nearest neighbor
of each of the $32\times 32$ crops, ordered in a grid representing the
location of the crop relative to the center crop.

(c) The four different training images retrieved as the nearest neighbors for
some crop of the image. The colored circles represent the corresponding color
in Figure 11(b).

(d) Four random images of a lens cap as a reference.

Figure 11: For every test image, we calculate DINOv2’s representation for each
of the $32\times 32=1024$ different crops. While (a) shows a monotonous
continuous translation of the cropped image in both dimensions, (b) shows that
the nearest neighbor alternates between several different images as shown in
(c). (d) serves as a reference, illustrating that while DINOv2’s
representation isn’t robust, it is meaningful in the sense that the nearest
neighbors are similar to the input image (a baby with the lens cap).

(a) Sixteen out of the possible $32\times 32=1024$ images of size $224\times
224$ pixels cropped from a $256\times 256$ image. The main object (a redshank
bird) is present in all the images.

(b) The index of the image in the training set, which is the nearest neighbor
of each of the $32\times 32$ crops, ordered in a grid representing the
location of the crop relative to the center crop.

(c) The four different training images retrieved as the nearest neighbors for
some crop of the image. The colored circles represent the corresponding color
in Figure 12(b).

(d) Four random images of a redshank as a reference.

Figure 12: For every test image, we calculate DINOv2’s representation for each
of the $32\times 32=1024$ different crops. While (a) shows a monotonous
continuous translation of the cropped image in both dimensions, (b) shows that
the nearest neighbor alternates between several different images as shown in
(c). (d) serves as a reference, illustrating that while DINOv2’s
representation isn’t robust, it is meaningful in the sense that the nearest
neighbors are similar to the input image (a redshank bird).

### Appendix 3 - More examples of DINOv2 ond Open-Clip change in predicted
nearest neighbor

(a)

(b)

(c)

(d)

Figure 13: More examples of the paper’s figure 2 - a change in the nearest
neighbor of both DINOv2 and Open-Clip as a function of a single pixel
translation. The main object is clearly seen in all of the test images

### Appendix 4 - More examples of AFC change in predicted nearest neighbor
with non-cyclic translations

(a)

(b)

(c)

### Appendix 5 - ImageNet-A Adversarial Robustness

Figure 15: When evaluated by predicted label, the adversarial robustness is
somewhat correlated the the accuracy, therefor, when evaluating on the hard
dataset of ImageNet-A (Hendrycks et al.) we can wee DINO is even less robust.
RICS improves the robustness, but it does decline as the translations becomes
bigger as the theory suggests.

### Appendix 6 - Robust Inferance by Crop Selection with DINOv2 results on the
PLACES365 Dataset

As stated in the paper, our method is proven to convert any classifier to a
robust classifier. We have shown it’s true for different models, and here we
show empirically it’s true for different dataset as well. We evaluate DINOv2
[23] with and without RICS on the PLACES365 [40] dataset which is more diverse
and challenging task than ImageNet [17].

Figure 16: Adversarial robustness (left) and Consistency (right) of the predicted nearest neighbor on the Places365 dataset. Using RICS improves the robustness dramatically in both metrics, for any shift size. Figure 17: Adversarial robustness of the predicted nearest neighbor for 1 pixel shift and model accuracy with KNN (k=11) on the Places365 dataset. Using RICS-rand improves the robustness dramatically with a 1% drop of accuracy only. Model | Accuracy |  Adv-Rob (1-NN) Shift Size: |  Consistency (1-NN) Shift Size:
---|---|---|---
| KNN | 1 | 3 | 5 | 9 | 1 | 3 | 5 | 9
DINOv2 | 49.03 | 74.94 | 54.77 | 45.98 | 37.04 | 88.83 | 83.15 | 80.15 | 77.20
RICS-MH (DINOv2) | 45.61 | 96.89 | 90.96 | 85.20 | 74.60 | 98.82 | 97.30 | 95.79 | 92.74
RICS-rand (DINOv2) | 47.97 | 96.68 | 90.53 | 84.90 | 74.55 | 98.75 | 97.15 | 95.65 | 92.76
Table 3: RICS is able to make DINOv2 much more robust, with a cost of 1% -
3.5% in accuracy. Results are calculated based on a random sample of 1000
images from the PLACES365 validation set.


# Adjusted composite likelihood
for robust Bayesian meta-analysis

Michele Lambardi di San Miniato  Department of Statistical Sciences,
University of Padova,
Via Cesare Battisti 241, 35121 Padova, Italy  Nicola Sartori  Department of
Statistical Sciences, University of Padova,
Via Cesare Battisti 241, 35121 Padova, Italy

(2021-03-20)

###### Abstract

A composite likelihood is a non-genuine likelihood function that allows to
make inference on limited aspects of a model, such as marginal or conditional
distributions. Composite likelihoods are not proper likelihoods and need
therefore calibration for their use in inference, from both a frequentist and
a Bayesian perspective. The maximizer to the composite likelihood can serve as
an estimator and its variance is assessed by means of a suitably defined
sandwich matrix. In the Bayesian setting, the composite likelihood can be
adjusted by means of magnitude and curvature methods. Magnitude methods imply
raising the likelihood to a constant, while curvature methods imply evaluating
the likelihood at a different point by translating, rescaling and rotating the
parameter vector. Some authors argue that curvature methods are more reliable
in general, but others proved that magnitude methods are sufficient to
recover, for instance, the null distribution of a test statistic. We propose a
simple calibration for the marginal posterior distribution of a scalar
parameter of interest which is invariant to monotonic and smooth
transformations. This can be enough for instance in medical statistics, where
a single scalar effect measure is often the target.

Keywords— inference; coverage; pseudo likelihood; clustered data.

##  1 Introduction

Likelihood functions are popular tools both in frequentist and Bayesian
inference. A genuine, or proper, likelihood is defined in terms of a
probability distribution, evaluated at the observed full dataset and treated
as a function of the parameter only. In the Bayesian framework, under a
suitable prior, a genuine likelihood implies a posterior distribution for the
parameter that can be used in inference. A composite likelihood is defined by
multiplying two or more genuine likelihood components together, but generally
this product is not a genuine likelihood in turn. This happens because the
product implies a working independence assumption that may not be valid in the
light of background knowledge. If each of the genuine likelihood components is
correctly specified, the maximum composite likelihood estimator is at least
consistent under a suitable asymptotic theory. Moreover, the estimation
variance can be assessed via a suitably defined sandwich covariance matrix.
For more on the theory of composite likelihoods in the classical framework,
one can refer to the review by Varin et al.,  (2011).

As to the Bayesian framework, a composite likelihood is not based on the
density of a correctly specified model, so it represents a wrong probability,
and the inverse probabilities implied by Bayes’ theorem are incorrect too:
this may imply either liberal or conservative inference, as first warned by
Monahan and Boos,  (1992), depending for instance on the sign of neglected
correlations among other things. For more on diagnostic methods to assess this
type of behavior, one can refer to Lazar,  (2003).

Consider the case of clustered data. In a rather simple case, the observations
will just make up data pairs. These pairs will be independent of each other
and internally correlated in general. This way of modeling can be useful in
clinical meta-analyses, where the treatment and control group are paired
within the studies and the studies are independent. A composite likelihood can
help to infer a marginal or conditional distribution, without the need to
assume an explicit form for the internal correlations. One can refer for
instance to Chandler and Bate,  (2007) on the topic of composite likelihoods
for clustered data.

Chen et al.,  (2015) consider the case of a meta-analysis on clustered count
data, where the counts of patients experiencing a given event are reported
among two treatment groups, say 1 and 2, for each study. Even if the groups
are independent by randomization, a spurious correlation between groups and
within studies can be observed as a result of heterogeneity between the
studies. See for instance Figure 1, where we report the logarithm of event
odds along with standard errors for each treatment group and study. We also
computed the summary log-odds, along with standard errors estimated in a
robust fashion accounting for correlations. The studies may be heterogeneous
just by involving patients with different baseline frailties: this is
reflected in the Figure by the presence of a covariant gradient among log-
odds.

Modeling the spurious correlation is not of substantial interest and it
essentially a side-effect of meta-analysis. Chen et al.,  (2015) resort to a
marginal composite likelihood, for which they specify simply the marginal
distribution for the outcome in each group. In order to assess estimation
variance they only have to assume the existence of the spurious correlation
within studies. They can only make inference on functions of marginal
parameters, which are actually their sole interest. This approach clearly
saves from making potentially wrong assumptions that are also difficult to
formulate, compare, and check in practice.

One can compare this approach with an alternative, also considered by Chen et
al.,  (2015), which involves explicit modeling of correlations within studies.
The most flexible approach in this sense is probably the copula modeling
approach, for which one may refer for instance to Nelsen,  (2007). It is
difficult though to find an appropriate copula model in every situation, as
there might also be too little studies for a meta-analysis, and this kind of
approach will likely require a model selection step. This is actually avoided
by means of a composite likelihood approach.

As to the Bayesian framework, a composite likelihood will usually require an
adjustment, or modification, that essentially preserves the composite’s
definition, while correctly accounting for uncertainty. Posterior
distributions could be mis-calibrated otherwise, as first warned by Monahan
and Boos,  (1992) for the more general case of pseudo likelihoods. Adjustments
are generally motivated on the basis of asymptotic theory, as in Miller,
(2019). The aim is to approximately recover some nice feature of genuine
posterior distributions, in the logic illustrated by Stoehr and Friel,
(2015). Ribatet et al.,  (2012) review two main classes of adjustments broadly
referred to as curvature and magnitude methods. As they make further
distinctions between overall and adaptive adjustments, we will only deal with
overall adjustments, which are simpler because they require less tuning. Set
aside their specific motivations, these methods can still be tested in
simulations, where the ground truth is known by design and several datasets
are simulated and analyzed automatically. Calibration plots can be used to
check the validity of probabilistic statements based on the adjusted posterior
distributions, in the fashion illustrated by Lazar,  (2003).

The most flexible adjustment is the curvature method, which involves rescaling
and rotating the parameter vector, in order to make posterior variance close
to the sandwich variance from the frequentist framework. The log-posterior is
asymptotically quadratic, so this is sufficient to recover approximate
calibration for all parameters. This method is advocated for by Ribatet et
al.,  (2012) when the analysis has no specific target but potentially aims at
making inference on the whole parameter vector.

Pauli et al.,  (2011) aim to recover the null chi-squared distribution of a
composite likelihood ratio test. In doing this, they first recall Bartlett
corrections, which essentially multiply the log-likelihood by a constant. This
constant is also called temperature and the adjustment itself is also known as
tempered composite likelihood. This operation results into a magnitude
adjustment, which can be enough to target specific aspects of the posterior
distribution, but it is still not flexible and it will not work under
reparameterization. The authors proposing it deem the method as targeted to
the whole parameter vector, so we will refer to it as an omnibus adjustment.

In this paper we argue that the omnibus magnitude adjustment used by Pauli et
al.,  (2011) is less safe than claimed, and we show this fact in the case of
an independence composite likelihood. Curvature methods advocated for by
Ribatet et al.,  (2012), for instance, constitute a safer alternative in
general, but in many applications the aim can be to just make inference on a
scalar parameter of interest: in such a case, we propose a simple targeted
magnitude adjustment that calibrates the marginal posterior distribution of
specifically that parameter of interest. Our adjustment is also invariant to
transformations of the target parameter.

We compare the omnibus magnitude adjustment with our targeted proposal and
with the curvature adjustment. As matrix square roots are involved in the
latter, we compare a definition of such square root that is common in
applications with an alternative definition proposed by Kessy et al.,  (2018):
this proposal reduces the implied rotation of the parameter vector implied by
curvature adjustments and preserves better some features of the composite
likelihood, including its asymmetries, skewness, kurtosis and so forth. A
simulation study on a clustered data model with an independence likelihood
will show that a magnitude adjustment can suffice in making inference about a
scalar parameter of interest: in this sense our adjustment seems as safe as a
curvature, though limited a priori to the target parameter, while Pauli et
al.,  (2011)’s version of magnitude adjustment implies no substantial
correction because of some mathematical subtleties that will be illustrated
later.

In medical statistics and other disciplines, where some effect measures
approximate each other, targeting a given parameter may also target
approximately a related one, making the adjustment reliable to a wider extent.
If there is a focus, targeted magnitude adjustments can be as good performers
as the curvature method, while remaining unreliable in general for all
parameters but the target.

##  2 Composite likelihood

Let $\theta$ be a $p$-dimensional parameter. Let $y=(y_{1},\dots,y_{K})$ be an
observed outcome of the data-generating process with true data distribution
$f_{0}(y)$. The outcome is modeled by a random vector $Y=(Y_{1},\dots,Y_{K})$.
We assume that one is willing to formulate
$P_{1}(\cdot;\theta),\dots,P_{K}(\cdot;\theta)$, the marginal distributions of
$Y_{1},\dots,Y_{K}$, respectively, whereas no interest lies on the joint
distributions. We assume that at least these marginals are correct, so there
exists a true parameter $\theta_{0}$ such that

$P_{k}(y_{k};\theta_{0})=\int f_{0}(y)\,dy_{-k}$

Then, $K$ genuine likelihoods $\mathcal{L}_{1},\dots,\mathcal{L}_{K}$ can be
defined, as follows:

$\mathcal{L}_{1}(\theta)=P_{1}(y_{1};\theta)\,,\quad\dots\,,\quad\mathcal{L}_{K}(\theta)=P_{K}(y_{K};\theta)\,.$

These definition can hold up to a proportionality constant. A composite
likelihood $\mathcal{L}_{C}(\theta)$ can be defined as their product

$\mathcal{L}_{C}(\theta)=\prod_{k=1}^{K}\mathcal{L}_{k}(\theta)\,.$

One may read $\mathcal{L}_{C}(\theta)$ in the expression above as a measure of
the probability to observe $y$ under $\theta$: then, the product in the
definition of $\mathcal{L}_{C}(\theta)$ implies the assumption that the random
vectors $Y_{1},\dots,Y_{K}$ are all independent. This is referred to as a
working independence assumption, which can be wrong without affecting the
consistency of estimators.

Let $\partial$ denote the gradient or Jacobian operator. We also need to
define the following likelihood-related quantities. Notation is mostly in line
with the review on composite likelihoods offered by Varin et al.,  (2011).

$\ell_{k}(\theta)=\log\mathcal{L}_{k}(\theta)$ is the generic $k$-th log-
likelihood.

$u_{k}(\theta)=\partial\ell_{k}(\theta)/\partial\theta$ is the generic $k$-th
score function.

$u(\theta)=\sum_{k=1}^{K}u_{k}(\theta)$ is the score function.

$\hat{H}_{k}(\theta)=-\partial u_{k}(\theta)$ is the generic $k$-th sample
sensitivity matrix.

$H_{k}(\theta)=\mathrm{E}\left\\{\hat{H}_{k}(\theta);\theta\right\\}$ is the
generic $k$-th sensitivity matrix.

$H(\theta)=\sum_{k=1}^{K}H_{k}(\theta)$ is the sensitivity matrix.

$J(\theta)=\mathrm{Var}\left\\{u(\theta);\theta\right\\}$ is the variability
matrix.

$\hat{J}(\theta)$ is a sample variability matrix that serves as a consistent
estimator of $J(\theta)$.

After a suitable definition of sample size $n$, under regularity conditions,
the maximum composite likelihood estimator $\hat{\theta}$ maximizing the
composite likelihood is a consistent estimator of $\theta$. The variance of
$\hat{\theta}$ can be approximated by the robust asymptotic variance
$\mathrm{V}_{\theta}$, defined as

$\mathrm{V}_{\theta}=H(\theta)^{-1}J(\theta)H(\theta)^{-1}\,.$

This is called robust, as contrasted to the naive variance
$\mathrm{N}_{\theta}$, which is defined as

$\mathrm{N}_{\theta}=H(\theta)^{-1}\,.$

The latter provides the variance assessment when trusting the composite
likelihood as if it were genuine. The robust and naive variances are defined
as the inverse of Godambe’s and Fisher’s information matrices, respectively.

In practice, $\mathrm{V}_{\theta}$ and $\mathrm{N}_{\theta}$ can be estimated
feasibly with

$\hat{\mathrm{V}}_{\theta}=\hat{H}(\hat{\theta})^{-1}\hat{J}(\hat{\theta})\hat{H}(\hat{\theta})^{-1}\quad\text{and}\quad\hat{\mathrm{N}}_{\theta}=\hat{H}(\hat{\theta})^{-1}\,.$

Assume a scalar parameter $\eta=\eta(\theta)$ is of interest, with
$\eta(\cdot)$ suitably smooth. Let $\partial\eta$ be the gradient of $\eta$
with respect to $\theta$. The maximum composite likelihood estimator is
equivariant, so $\hat{\eta}=\eta(\hat{\theta})$. Moreover, $\hat{\theta}$ is
consistent and has bias of order $O(n^{-1})$, and the same is true for
$\hat{\eta}$. The robust variance of $\hat{\eta}$, denoted by
$\mathrm{V}_{\eta}$, can be approximated using the delta method, with

$\mathrm{V}_{\eta}=\partial\eta^{\top}\cdot\mathrm{V}_{\theta}\cdot\partial\eta\,.$

We call this the robust delta variance of $\hat{\eta}$. A naive version can be
defined as

$\mathrm{N}_{\eta}=\partial\eta^{\top}\cdot\mathrm{N}_{\theta}\cdot\partial\eta\,.$

Estimating the naive and the robust variance for $\eta$ also involves
evaluating $\partial\eta$ at $\hat{\eta}$.

###  2.1 Clustered data case

The definition of $\hat{J}(\theta)$ can be based on background knowledge about
the true correlation structure, which is neglected under the working
independence assumption when defining the composite likelihood. Consider a
meta-analysis summarizing $N$ studies: each study implies a cluster of
observations; the clusters are independent of each other, but they are to be
treated as internally correlated so to account for heterogeneity-related
issues.

One can define $Y_{k}=(Y_{1k},\dots,Y_{ik},\dots,Y_{Nk})$, for $k=1,\dots,K$,
where the studies are indexed by $i=1,\dots,N$ and the random variable
$Y_{ik}$ models the outcome $y_{ik}$ from the $k$-th treatment group of the
$i$-th study.

The following likelihood-related quantities can be defined.

$P_{ik}(\cdot;\theta)$ is the distribution of $Y_{ik}$,

$\mathcal{L}_{ik}(\theta)=P_{ik}(y_{ik};\theta)$ is the $i$-th study’s
contribution to $\mathcal{L}_{k}(\theta)$, so
$\mathcal{L}_{k}(\theta)=\prod_{i=1}^{N}\mathcal{L}_{ik}(\theta)$,

$\ell_{ik}(\theta)=\log\mathcal{L}_{ik}(\theta)$ is the $i$-th study’s
contribution to $\ell_{k}(\theta)$, so
$\mathcal{L}_{k}(\theta)=\sum_{i=1}^{N}\ell_{ik}(\theta)$,

$u_{ik}(\theta)=\partial\ell_{ik}(\theta)$ is the $i$-th study’s contribution
to $u_{k}(\theta)$, so $u_{k}(\theta)=\sum_{i=1}^{N}u_{ik}(\theta)$,

$u_{i*}(\theta)=\sum_{k=1}^{K}u_{ik}(\theta)$ is the $i$-th study’s
contribution to $u(\theta)$, so $u(\theta)=\sum_{i=1}^{N}u_{i*}(\theta)$.

Just by assuming the correlation structure, the following holds:

$J(\theta)=\mathrm{Var}\left\\{u(\theta);\theta\right\\}=\mathrm{Var}\left\\{\sum_{i=1}^{N}u_{i*}(\theta);\theta\right\\}=\sum_{i=1}^{N}\mathrm{Var}\left\\{u_{i*}(\theta);\theta\right\\}\,.$

So, by treating $u_{i*}(\theta)$ as a column vector, an estimator for
$J(\theta)$ can be defined as

$\hat{J}(\theta)=\sum_{i=1}^{N}u_{i*}(\theta)\cdot u_{i*}(\theta)^{\top}\,.$

##  3 Bayesian adjustment

The composite likelihood approach has clear advantages over full likelihoods
in the case we consider because there could hardly be any knowledge about the
appropriate copula model for the problem at hand. The composite likelihood
does not need to explicitly formulate such uninteresting aspects of the model.

Composite likelihoods are more naturally dealt with in the classical framework
than the Bayesian one. Monahan and Boos,  (1992) discussed the risks of mis-
calibrated posterior distributions when dealing with pseudo, not genuine,
likelihoods, such as the pseudo likelihood from Cox’ proportional hazards
model. A composite likelihood may be replaced with an artificial likelihood in
the Bayesian framework, as proposed by Müller,  (2013): the artificial
likelihood is Gaussian, with mean vector equal to the maximum composite
likelihood estimator $\hat{\theta}$ and covariance matrix equal to the robust
variance $\mathrm{V}_{\theta}$. This artificial likelihood is also called the
sandwich likelihood and it can also be motivated under the asymptotic theory
developed by Miller,  (2019).

A composite likelihood is not expected to be well calibrated and therefore
needs to be adjusted, or slightly modified, before using it in the Bayesian
framework. The adjustment’s tuning can be based on the output of a classical
analysis and the resulting adjusted composite likelihood should imply an
asymptotically calibrated posterior distribution. A reasonable adjustment
method should be defined so that one approximately falls back to the composite
likelihood when the working independence assumption is correct. Another
requirement is that the adjustment should preserve, as far as possible, the
shape and the asymmetries of the original composite likelihood: the sandwich
likelihood by Müller,  (2013) is not particularly suited on this aspect.

The main types of adjustments discussed by Ribatet et al.,  (2012) are
curvature and magnitude methods. The curvature method transforms the parameter
vector before feeding it to the composite likelihood, while the magnitude
transforms the composite likelihood itself. The magnitude transformation
raises the likelihood to a constant or, equivalently, multiplies the log-
likelihood by the same constant. The curvature transformation is affine or
linear. Actually, the magnitude method involves just one tuning constant, as
contrasted to the curvature method, which involves a tuning matrix made up of
$O(p^{2})$ constants. This is reflected in the higher flexibility of curvature
method with respect to magnitude methods.

###  3.1 Curvature adjustment

The curvature method replaces the composite likelihood $\mathcal{L}_{C}$ with
the curvature-adjusted composite likelihood ${\mathcal{L}_{C}}_{adj}$, defined
as

${\mathcal{L}_{C}}_{adj}(\theta;A)=\mathcal{L}_{C}(A(\theta-\hat{\theta})+\hat{\theta})\,.$

In the expression above, $A$ is a $p\times p$ tuning matrix. In the terms
proposed by Stoehr and Friel,  (2015), it can be seen as a device that allows
the adjusted composite likelihood to satisfy a matching rule stated as
follows:

$\mathrm{N}_{\theta}({\mathcal{L}_{C}}_{adj}(\cdot;A))=\mathrm{V}_{\theta}\,.$

Here, $\mathrm{N}_{\theta}({\mathcal{L}_{C}}_{adj})$ is the naive variance
matrix for $\theta$, as implied by the adjusted composite likelihood
${\mathcal{L}_{C}}_{adj}$. This aspect of ${\mathcal{L}_{C}}_{adj}$ and the
extreme point $\hat{\theta}$ are the only features of
${\mathcal{L}_{C}}_{adj}$ relevant to the asymptotic theory as developed, for
instance, by Miller,  (2019). The matrix $A$ that feasibly solves the problem
above is in the form:

$A=\left\\{(\hat{\mathrm{N}}_{\theta}^{-1})^{1/2}\right\\}^{-1}(\hat{\mathrm{V}}_{\theta}^{-1})^{1/2}\,.$

Here, $B=A^{1/2}$ denotes the square root of matrix $A$, such that
$B^{\top}B=A$. As a matter of linear algebra, the matrix square root is not
unique: for every orthonormal matrix $E$, also $E\cdot B$ can serve as a
square root of $A$. This phenomenon is called rotational freedom. Different
definitions of the square root imply different optimality or invariance
properties, as discussed by Kessy et al.,  (2018).

In simulations we considered a main version of curvature adjustment based on a
definition of matrix square root referred to as zero-phase component analysis
(ZCA) in Kessy et al.,  (2018). This definition is popular among
practitioners, even under a different name, since it has the same eigen-
decomposition as the argument matrix excepted the eigenvalues, which are
square-rooted. ZCA is also symmetric whenever the argument is, and it is
endowed with the nice property

$(A^{1/2})^{-1}=(A^{-1})^{1/2}\,.$

We also considered an experimental variant of curvature adjustment, involving
the ZCA-cor definition of matrix square root, which is invariant to rescaling
operations when the argument is a covariance matrix. This alternative implied
no substantial improvement in simulations but it is reported for reference.
This should be a concern in general, though, because different matrix square
roots imply different rotations that can more or less alter the asymmetries of
the composite likelihood. Such a rotation effect does not matter
asymptotically, but it can do in finite samples.

Shaby,  (2014) proposed an adjustment similar to curvature method in some
respect, but it first draws samples from the unadjusted posterior distribution
and then scales and rotates these samples around $\hat{\theta}$ so that its
covariances will reflect the robust variance $\mathrm{V}_{\theta}$. This
solution meets the asymptotic requirements for Miller,  (2019) but it reduces
the contribution of the prior to the analysis. It will not be discussed
further in this paper, though it could be a useful alternative to curvature
adjustments from an objective Bayesian perspective.

###  3.2 Magnitude adjustments

The curvature-adjusted posterior is approximately calibrated for all
parameters and thus also for any target $\eta$, but this solution is
excessively flexible if just $\eta$ is of interest. One may consider an
adjustment that calibrates for $\eta$ without calibrating for other
parameters, especially if this approach requires less constants to be tuned,
resulting into improved stable behavior.

One can magnitude-adjust the composite likelihood by raising it to a suitable
power $w$, called temperature, in the following fashion.

${\mathcal{L}_{C}}_{adj}(\theta;w)=\mathcal{L}_{C}(\theta)^{w}\,.$

The term magnitude refers to the change in the order of magnitude of the
relative log-likelihood at any given $\theta$ implied by $w$.

This formulation arises naturally for instance in the work by Pauli et al.,
(2011), as a Bartlett correction, as their goal is to recover the null
distribution of some likelihood ratio test. In particular, they set

$w=\frac{p}{\mathrm{trace}\left\\{(\mathrm{N}_{\theta})^{-1}\mathrm{V}_{\theta}\right\\}}\,,$

where the unary $\mathrm{trace}$ operator defined for square-sized matrices
returns the sum of the argument’s diagonal entries. This is referred to as the
omnibus magnitude adjustment in the following. In general, even with a
different aim, the tuning constant $w$ allows for less corrections than
implied by curvature adjustment’s tuning matrix $A$, but this is actually not
a problem if the goal is suitably limited by some specific inferential needs.

Under the same logic seen before, as outlined by Stoehr and Friel,  (2015),
consider a magnitude-adjusted composite likelihood ${\mathcal{L}_{C}}_{adj}$
that satisfies a different matching rule, which is less demanding and stated
as follows:

$\mathrm{N}_{\eta}({\mathcal{L}_{C}}_{adj}(\cdot;w))=\mathrm{V}_{\eta}\,.$

Here, $\mathrm{N}_{\eta}({\mathcal{L}_{C}}_{adj})$ is the naive variance for
$\eta$ alone, as implied by the adjusted composite likelihood
${\mathcal{L}_{C}}_{adj}$. The tuning constant $w$ that feasibly solves the
problem above is:

$w=\frac{\mathrm{N}_{\eta}}{\mathrm{V}_{\eta}}\,.$

The composite log-likelihood is thus rescaled in order to imply the robust
variance at least for the target parameter $\eta$. The posterior will be
approximately calibrated for $\eta$ after Miller,  (2019). The proposed
magnitude adjustment is based on delta variances, so it is the same by design
for all smooth invertible transformations of $\eta$.

##  4 Marginal beta-binomial likelihood

Chen et al.,  (2015) consider the case of a meta-analysis with independence
composite likelihood on count data: the number of patients experiencing a
given event of interest is reported for each treatment group in each study.
The most natural choice is to model the counts as binomial-distributed, but
over-dispersion is often observed in practice, as a result of heterogeneity
among studies.

It is convenient to resort to hierarchical modeling. Conditional to the
probability of event $p_{ik}$ in the $k$-th group of the $i$-th study,
$Y_{ik}$ is binomial-distributed with probability of event $p_{ik}$ and known
size $n_{ik}$, and it is independent of all other observable outcomes in the
light of randomization. If there were no heterogeneity, it would hold
$p_{ik}=\pi_{k}$ for all $i=1,\dots,N$, but the probabilities can otherwise be
assumed random, with expected value $\pi_{k}$ independent of the study.

As detailed by Chen et al.,  (2015), a marginal Beta distribution with shape
parameters $\alpha_{k},\beta_{k}$ can be assumed for $p_{ik}$, so the average
probability of event in the $k$-th treatment group be defined as

$\mathrm{E}(p_{ik};\theta)=\pi_{k}=\frac{\alpha_{k}}{\alpha_{k}+\beta_{k}}\,,$

where $\theta=(\alpha_{1},\beta_{1},\dots,\alpha_{K},\beta_{K})$. Beta models
can flexibly represent random quantities bounded between $0$ and $1$. This
choice has some other computational advantages since it implies that,
marginally, thus without conditioning on $p_{ik}$, $Y_{ik}$ is beta-binomial-
distributed. In formulas, it holds

$P_{ik}(y;\theta)=p(y;n_{ik},\alpha_{k},\beta_{k})=\binom{n_{ik}}{y}\frac{B(y+\alpha_{k},n_{ik}-y+\beta_{k})}{B(\alpha_{k},\beta_{k})}\,,$

where $B(\cdot,\cdot)$ is the beta function.

Conditional on probabilities, the counts can be assumed independent due to
randomization, but the observed correlations within studies are accounted for
in terms of a correlation among latent probabilities $p_{i1},\dots,p_{iK}$.
This correlation does not need to be specified for the composite likelihood
approach to be used. Actually, the approach is useful in that it does not
require to specify that aspect of the model. One can imagine though that there
exists a copula model that rules the probabilities identically within each
study, resulting in beta marginals for the probabilities. After Sklar’s
theorem, as per Nelsen,  (2007), such a copula exists and is unique, because
the beta distribution is absolutely continuous.

An average odd of event in the $k$-th treatment group can be defined as

$o_{k}=\frac{p_{k}}{1-p_{k}}\,.$

Chen et al.,  (2015) consider the case where $K=2$ and they aim at making
inference on a scalar parameter of interest called diagnostic odds ratio
($\rho$), which is a function of the marginal parameters defined as:

$\rho=\frac{o_{1}}{o_{2}}\,.$

This quantity generalizes the odds ratio when probabilities are random.

In our study, we also consider an alternative effect measure, defined as

$\delta=p_{1}-p_{2}\,.$

The two measures $\rho$ and $\delta$ are related and legitimate alternatives
as an effect measure of interest. Our concern is whether an adjusted composite
likelihood that is targeted at $\rho$ is also reliable for $\delta$. This may
not be the case in switching to the Bayesian framework, where the composite
likelihood can be magnitude-adjusted in order to recover calibration for
potentially few aspects of the model.

##  5 Calibration plots

Beyond the rationales, different adjustments can be compared in terms of
approximate recovery of posterior calibration: this is the property of a
posterior distribution that yields asymptotically valid probabilistic
statements. In repeated simulations for which the ground truth is known, the
posterior quantile of order $p$ for a given parameter can be compared with the
true value of the parameter, which should lie under the posterior quantile
with frequency about $p$. The importance of these concerns will naturally
depend on the way of interpreting probability. From a classical perspective,
Miller,  (2019) predicts an asymptotic agreement between classical and
Bayesian analyses, so a total disagreement between the two might be
concerning.

Lazar,  (2003) dealt with empirical likelihoods and assessed posterior
calibration for a single scalar parameter of interest. The same approach can
be applied to composite likelihood functions, as the aim is the same in
checking for posterior calibration. In simulations, the ground truth
$\theta_{0}$ and so $\eta_{0}=\eta(\theta_{0})$ are known, and some properties
of the posterior $\mathrm{Pr}(\cdot\mid Y=y)$ can be assessed. To this end,
Lazar,  (2003) defined a statistic, which we will denote $h$. The statistic
$h$ for $\eta$ is defined as

$h=\mathrm{Pr}(\eta\leq\eta_{0}\mid Y=y)\,.$

Under repeated sampling, asymptotically it holds

$h\,\dot{\sim}\,\mathcal{U}\mathrm{niform}[0,1]\,.$

We define the statistic $g$ as

$g=|2h-1|\,,$

which has the same uniform limiting distribution. We use quantile-based
credible intervals (QB) instead of highest-posterior-density intervals (HPD)
for the sake of simplicity and of invariance to monotonic transformations of
$\eta$. If the credibility level is set equal to $p$, the QB interval for
$\eta$ is delimited by $\mathrm{q}(\eta;p/2\mid Y)$ and
$\mathrm{q}(\eta;1-p/2\mid Y)$. Here, $\mathrm{q}(\eta;p\mid Y)$ is $\eta$’s
posterior quantile of order $p$, such that

$\mathrm{Pr}\left\\{\eta\leq\mathrm{q}(\eta;p\mid Y)\mid Y\right\\}=p\,.$

In this paper we call calibration plot the empirical cumulative distribution
function of statistic $g$. This chart can be read in terms of the effective
coverage, along the $y$-axis, as a function of the nominal coverage, along the
$x$-axis, of QB credible intervals. The ideal behavior for the curve is a
straight line connecting points $(0,0)$ and $(1,1)$, which would imply perfect
calibration. Curves above this line denote conservative inference, while any
curve below would denote liberal inference.

One must consider that the ideal behavior may not be achieved in finite
samples, due to non-Gaussian likelihood or informative priors to begin with,
so the performance of curvature methods will set the standard in the light of
the arguments made by Ribatet et al.,  (2012).

##  6 Simulation study

In simulations, we generate 1000 datasets for each of some selected settings.
Every setting is a joint configuration of copula model, its implied
correlation, and a phase, all defined as follows.

The copula model, rules the latent probabilities within studies, and ranges in
Clayton, Frank, Gumbel copulas. The copula model implies a rank correlation
among probabilities within studies, ranging from $0.5$ to about $0.95$ in our
simulations; it is set by suitably tuning the copula’s only parameter.

Some popular parametric copula models were chosen. These are available in the
R package copula. Under these models, the family’s only parameter determines
rank correlations in a identifiable way, independent of marginal
distributions. This relation is also exploited by Tsukahara,  (2005) for semi-
parametric estimation of copulas, because it does not require the estimation
of marginal distributions.

The marginal parameters are

$\theta=(\log\alpha_{1},\log\beta_{1},\log\alpha_{2},\log\beta_{2})\,,$

so the gradients of $\rho$ and $\delta$ to compute delta variances are

$\frac{\partial\rho}{\partial\theta}=\frac{1}{\rho}\cdot(1,-1,-1,1)\quad\text{and}\quad\frac{\partial\delta}{\partial\theta}=\left[\pi_{1}\cdot(1,-1),\pi_{2}\cdot(-1,1)\right]\,.$

The logarithm is used in the definition of $\theta$ as the posterior is
supposed to be better approximated by a Gaussian when the parameters can take
any real value. This choice is useful in the Metropolis-Hastings sampling of
the posterior distribution, of which we use an adaptive variant implemented in
the R package MHadaptive. The same parameterization is also used in the R
package xmeta for the assessment of estimation variance. Burn in and thinning
were used to reduce chains’ autocorrelation and improve on mixing, which made
the simulations computationally burdening.

The sample sizes $n_{ik}$ for each treatment group and study are the same from
the meta-analysis dataset used by Chen et al.,  (2015), also available is R
package xmeta, while the data $y_{ik}$ were simulated.

As to phases, in phase 1 each dataset is generated conditional to marginal
parameters, which are randomly drawn for each dataset from an informative
prior for $\theta$. The prior chosen is an uncorrelated Gaussian distribution
with zero means and unit variances. In phase 2, the true marginal parameters
are fixed and they can take, equally likely, any of the two following
configurations:

$\exp(\theta^{\prime})=(3.11,2.91,3.94,3.36)\,,\quad\exp(\theta^{\prime\prime})=(0.5,0.5,0.5,0.5)\,.$

These configurations are the same used in experiments by Chen et al.,  (2015).
The prior in phase 2 is only used in the Bayesian analysis, and not to draw
$\theta$, and it is the same prior as in phase 1 but made uninformative by
setting all marginal variances equal to $10^{4}$ instead of $1$.

The configuration $\theta^{\prime}$ is chosen on purpose to be close to the
actual estimates obtained on the original data by Chen et al.,  (2015). The
concern is whether one can consistently infer the truth at least in a
neighborhood of some relevant parameter configuration, since the parameter
space cannot be probed exhaustively.

###  6.1 Results

At most 1 every 1000 simulations was troubled for either the classical
estimation step or the subsequent Bayesian analysis and was discarded. The
results reported in the following are thus conditional to no error reported in
code. Results are shown in Figures 2, 3, 4.

The unadjusted posterior here is equivalent to a magnitude adjustment with
$w=1$ and it is compared with the standard according to Ribatet et al.,
(2012), which is the curvature method. The standard is presented in two
variants, based on the definitions ZCA and ZCA-cor of matrix square roots
illustrated by Kessy et al.,  (2018). The magnitude-adjusted posterior is used
in both the omnibus version by Pauli et al.,  (2011) and our variant targeted
at $\rho$. Performance on a second parameter $\delta$ of potential interest is
monitored as well.

Each distinct curve in calibration plots summarizes either 1000 analyses from
phase 1 or 500 analyses from phase 2, the latter being drawn from
$\theta^{\prime}$ or $\theta^{\prime\prime}$, while the former cannot be
separated according to the ground truth, which is distinct each time by
sampling. The curves are colored according to the copula’s rank correlation,
since this seems to explain most of the variability.

As one can see, there are visible issues in three cases: when the unadjusted
posterior or the omnibus magnitude adjustment are used to infer $\eta$, and
when the magnitude adjustment targeted to $\rho$ is used to infer a marginal
parameter.

As to the first issue about the unadjusted posterior, when making inference on
a function $\eta$ of either $(\alpha_{1},\beta_{1})$ or
$(\alpha_{2},\beta_{2})$ alone, the relevant likelihood is just one of the
proper likelihood components, so the unadjusted composite likelihood yields
valid inference on $\eta$, but if $\eta$ is a function that involves both
groups of parameter components the information is assessed incorrectly.

As to the second issue about the omnibus magnitude adjustment, one has to
consider that here each marginal distribution is correctly specified and ruled
by a different set of parameters, so the matrix $H$ is block-diagonal and its
diagonal blocks coincide with the corresponding entries in $J$. So, after a
little linear algebra, one can see that the omnibus tuning constant is exactly
$w=1$ when $J$ is known. The equality holds approximately when $J$ is
estimated. This is especially concerning, because the adjustment does not work
even under the model.

As to the last issue about the targeted magnitude adjustment, the arguments by
Ribatet et al.,  (2012) hold, as they stress the unreliability of magnitude
adjustments as to anything but the target, see Figure 2. Here we stress the
fact that magnitude adjustments seem able to recover posterior calibration to
essentially the same extent as the curvature method at least as far as the
target $\rho$ is concerned, which is clearly the focus in the problem at hand,
see Figure 3.

As a last remark, parameters that are related to the target $\rho$, such as
$\delta$, can benefit from magnitude adjustments even when they are not
directly targeted, see Figure 4. It is useful to consider the expression for
gradients of $\rho$ and $\delta$: when the average probabilities $\pi_{1}$ and
$\pi_{2}$ are close to each other, those gradients are proportional and the
magnitude adjustment targeted at $\rho$ and the one targeted at $\delta$ are
similar, so one can target either of the two parameters of interest directly,
and this will approximately imply targeting the other too.

##  7 Discussion

In this paper we considered some Bayesian adjustment methods for composite
likelihoods as reviewed by Ribatet et al.,  (2012). They proposed an adaptive
version to these adjustments, along with the version we used for our study:
our version is called overall adjustment in their work and it is tuned once
before sampling the posterior. Their adaptive version essentially repeats the
tuning of constants for each method before each iteration in posterior
sampling à-la Gibbs. We did not use this adaptive version of the adjustments
as checking for convergence would have been more difficult. Convergence
diagnostics are not reported here, but they were checked in just few cases in
order to tune the burn in and thinning behavior of the posterior sampler.

The adjustments discussed here have some intuitive properties that make their
tuning rather straightforward. These methods require to perform a frequentist
analysis before the Bayesian one. Magnitude adjustments additionally require
stating a scalar target parameter, which a priori is the only parameter for
which the posterior will be approximately calibrated.

We hinted at how the analysis of clustered data with independence likelihoods
can benefit from magnitude adjustments in targeting specific parameters of
interest. Actually, we managed to make reliable inference even on a parameter
that was not directly targeted but still related and potentially relevant.
Some known effect measures in medical statistics can approximate each other,
so a magnitude adjustment can be useful to a broader extent than the nominal.
The additional flexibility of curvature adjustments is not needed in such
cases.

The warnings expressed by Ribatet et al.,  (2012) still stand, as magnitude
adjustments do not recover full posterior calibration, both in the omnibus and
the targeted version, but more in general under any tuning rule due to
intrinsic lack of flexibility. When the scope of the analysis is not limited
enough, the more flexible curvature adjustment should be used instead of the
omnibus magnitude method, which is generally unsafe.

## References

  * Chandler and Bate, (2007) Chandler, R. E. and Bate, S. (2007).  Inference for clustered data using the independence loglikelihood.  Biometrika, 94:167–183.
  * Chen et al., (2015) Chen, Y., Hong, C., Ning, Y., and Su, X. (2015).  Meta-analysis of studies with bivariate binary outcomes: a marginal beta-binomial model approach.  Statistics in Medicine, 35:21–40.
  * Chivers, (2012) Chivers, C. (2012).  MHadaptive: General Markov Chain Monte Carlo for Bayesian inference using adaptive Metropolis-Hastings sampling.  R package version 1.1-8.
  * Hofert et al., (2020) Hofert, M., Kojadinovic, I., Maechler, M., and Yan, J. (2020).  copula: Multivariate dependence with copulas.  R package version 1.0-1.
  * Hong et al., (2021) Hong, C., Luo, C., Tong, J., and Chen, Y. (2021).  xmeta: A toolbox for multivariate meta-analysis.  R package version 1.3-0.
  * Kessy et al., (2018) Kessy, A., Lewin, A., and Strimmer, K. (2018).  Optimal whitening and decorrelation.  The American Statistician, 72:309–314.
  * Lazar, (2003) Lazar, N. A. (2003).  Bayesian empirical likelihood.  Biometrika, 90:319–326.
  * Miller, (2019) Miller, J. W. (2019).  Asymptotic normality, concentration, and coverage of generalized posteriors.  arXiv preprint, arXiv:1907.09611v1.
  * Monahan and Boos, (1992) Monahan, J. F. and Boos, D. D. (1992).  Proper likelihoods for bayesian analysis.  Biometrika, 79:271–278.
  * Müller, (2013) Müller, U. K. (2013).  Risk of Bayesian inference in misspecified models, and the sandwich covariance matrix.  Econometrica, 81:1805–1849.
  * Nelsen, (2007) Nelsen, R. B. (2007).  An Introduction to Copulas.  Springer.
  * Pauli et al., (2011) Pauli, F., Racugno, W., and Ventura, L. (2011).  Bayesian composite marginal likelihoods.  Statistica Sinica, 21:149–164.
  * Ribatet et al., (2012) Ribatet, M., Cooley, D., and Davison, A. C. (2012).  Bayesian inference from composite likelihoods, with an application to spatial extremes.  Statistica Sinica, 22:813–845.
  * Shaby, (2014) Shaby, B. A. (2014).  The open-faced sandwich adjustment for MCMC using estimating functions.  Journal of Computational and Graphical Statistics, 23:853–876.
  * Stoehr and Friel, (2015) Stoehr, J. and Friel, N. (2015).  Calibration of conditional composite likelihood for bayesian inference on gibbs random fields.  arXiv preprint, arXiv:1502.01997v2.
  * Tsukahara, (2005) Tsukahara, H. (2005).  Semiparametric estimation in copula models.  Canadian Journal of Statistics, 33:357–375.
  * Varin et al., (2011) Varin, C., Reid, N., and Firth, D. (2011).  An overview of composite likelihood methods.  Statistica Sinica, 21:5–42.

##  Appendix A Plots

###  A.1 Dataset

A presentation on the dataset taken from Chen et al.,  (2015). Each axis
represents one out of two treatment groups, the logarithm of event odd in each
group is reported for each study in their meta-analysis. The summary for each
group has been estimated based on their marginal beta-binomial approach and
the standard error is based on the sandwich covariance matrix.

Figure 1:

###  A.2 Calibration plots

In the following pages, calibration plots are listed. For each combination of
phase, true $\theta$ if phase is 2, true copula and true correlation within
study, the empirical cumulative distribution function of statistic $h$ is
reported as a distinct line. Each line pools 1000 simulations if phase is 1
and 500 simulations if phase is 2. This is repeated for each adjustment method
considered in the study. The plot can be read as effective credibility level
(vertically) versus nominal credibility level (horizontally).

Lines are colored according to true correlation, which explains most of the
variability when this is present. It is evident that the higher the
correlation within studies, the worse the coverage properties of the
unadjusted posterior when looking at parameters $\rho$ and $\delta$, which are
functions of all components of $\theta$. In the trivial case where
probabilities are independent within studies, the composite likelihood is
proper and one may expect all methods to fall back to it.

Figure 2:

Figure 3:

Figure 4: