#  _DR 2S_: Deep Regression with Region Selection for Camera Quality
Evaluation

Marcelin Tworski  LTCI, Télécom Paris,
Institut polytechnique de Paris
Palaiseau, France
and DXOMARK
Boulogne-Billancourt, France
<EMAIL_ADDRESS>Stéphane Lathuilière  LTCI, Télécom Paris,
Institut polytechnique de Paris
Palaiseau, France
<EMAIL_ADDRESS>Salim Belkarfa  DXOMARK
Boulogne-Billancourt, France
<EMAIL_ADDRESS>Attilio Fiandrotti  LTCI, Télécom Paris,
Institut polytechnique de Paris
Palaiseau, France
<EMAIL_ADDRESS>Marco Cagnazzo  LTCI, Télécom Paris,
Institut polytechnique de Paris
Palaiseau, France
<EMAIL_ADDRESS>
###### Abstract

In this work, we tackle the problem of estimating a camera capability to
preserve fine texture details at a given lighting condition. Importantly, our
texture preservation measurement should coincide with human perception.
Consequently, we formulate our problem as a regression one and we introduce a
deep convolutional network to estimate texture quality score. At training
time, we use ground-truth quality scores provided by expert human annotators
in order to obtain a subjective quality measure. In addition, we propose a
region selection method to identify the image regions that are better suited
at measuring perceptual quality. Finally, our experimental evaluation shows
that our learning-based approach outperforms existing methods and that our
region selection algorithm consistently improves the quality estimation.

##  I Introduction

With the rapid rise of smartphone photography, many people got more interested
in their photographs and in capturing high-quality images. Smartphone makers
responded by putting an increasing emphasis on improving their cameras and
image processing systems. In particular, sophisticated imaging pipelines are
now used in the Image Signal Processor (ISP) of any smartphone. Such software
makes use of recent advances such as multi-camera fusion [1], HDR+ [2] or deep
learning-based processing to improve image quality. Despite the high
performance of those techniques, their increasing complexity leads to a
correspondingly increased need for quality measurement.

Several attributes are important to evaluate an image: target exposure and
dynamic range, color (saturation and white balance) texture, noise and various
artifacts that can affect the quality of the final image [3]. When it comes to
camera evaluation, a camera can be evaluated for its capabilities in low-
light, zoom, or shallow depth-of-field simulation.

In this work, we are interested in estimating the quality of a camera for
multiple lighting conditions. More specifically, we aim at evaluating camera
capabilities to preserve fine texture details. We also refer to this problem
as texture quality estimation. Texture details preservation should be measured
in a way that reflects human perception. A typical way to evaluate the quality
of a set of cameras consists of comparing shots of the same visual content in
a controlled environment. The common visual content is usually referred to as
a _chart_. The motivation for using the same chart when comparing different
cameras is twofold. First, it facilitates the direct comparison of different
cameras. In particular, when it comes to subjective evaluation, humans can
more easily provide pairwise preferences than an absolute quality score.
Second, when the common noise-free chart is known, this reference can be
explicitly included in the quality measurement process. In this context, the
Modulation Transfer Function (MTF) is a widely-used tool for evaluating the
perceived sharpness and resolution, which are essential dimensions of texture
quality. First, MTF-based methods suffer from important drawbacks. MTF-based
methods are originally designed for conventional optic systems that can be
modeled as linear. Consequently, non-linear processing in the camera
processing pipeline, such as multi-images fusion or deep learning-based image
enhancement, may lead to inaccurate quality evaluation [4].

Second, these methods assume that the norm of the device transfer function is
a reliable measure of texture quality. However, in this work, we argue that
this assumption fails to account for many nuances of human perception. Some
recent works have shown that the magnitude of image transformations do not
always coincide with the perceived impact of the transformations [5].
Consequently, we advocate that human judgment should more explicitly be
included in the texture quality measurement process.

As a consequence, we propose _DR 2S_, a Deep Regression method with Region
Selection. Our contributions are threefold. First, we formulate the problem of
assessing the texture quality of a device as a regression problem and we
propose to use a deep convolutional network (ConvNet) for estimating quality
scores. We aim at obtaining a score which would be close to a subjective
quality judgment: to this end, we use annotations provided by expert human
observers as ground-truth at training time. Second, we propose an algorithm to
identify the regions in a chart that are better suited to measure perceptual
quality. Finally, we perform an extensive evaluation study that shows that our
learning-based approach performs better than existing methods for texture
quality evaluation and that our region selection algorithm leads to further
improvement in texture quality measurement.

##  II Related Works

In this section, we review existing works on texture quality assessment.
Existing methods can be classified into two main categories: MTF-based and
learning-based methods.

###  II-A MTF-based methods

A simple model for a camera consists in a linear system that produces an image
$y$ as a convolution of the point spread function $h$ and the incoming radiant
flux $x$. In the frequency domain, $Y(f)=H(f)X(f)+N(f)$, where we also
consider additive noise $N$. The modulation transfer function is
$\operatorname{MTF}(f)=|H(f)|$ and it is commonly used to characterize an
optic acquisition device [6].

Acutance is a single-value metric calculated from the device MTF. To compute
this value, a contrast sensitivity function (CSF), modeling the spatial
response of the human visual system, is used to weight the values of the MTF
for the different spatial frequencies. The CSF depends on a set of parameters
named viewing conditions, which describes how an observer would look at the
resulting image. These parameters are usually the printing height and the
viewing distance. The acutance is defined as the ratio of the integral of the
weighted MTF by CSF’s integral.

The key assumption in MTF-based methods is that they dispose of the noise-free
content to estimate the transfer function. The noise-free content is often
referred to as _full reference_. Therefore, these methods are usually used
only with synthetic visual charts. To the best of our knowledge, the only work
estimating acutance on natural scenes was proposed by Van Zwanenberg et al.
[4]. This method uses edges detection in the scene to compute the MTF. Early
methods use charts containing a blur spot or a slanted edge for this
computation. Loebich et al. [7] propose a method using the Siemens-Star. Cao
et al. propose to use the Dead-Leaves model [8, 9], and introduce an
associated method in [10], which is shown to be more appropriate to describe
fine detail rendering since its texture are more challenging for devices. This
chart is usually referred to as the _Dead-Leaves_ chart. In this chart, the
reference image is generated with occluding disks generated with a random
center location, radius and grey-scale value. Importantly, digital camera
systems present high-frequency noise, which affects the MTF estimation by
dominating signal in the higher frequencies. Consequently, estimating the
noise power spectral density (PSD) is key to obtain an accurate acutance
evaluation, and this task is not easily performed on the textured region. As a
consequence, noise PSD is typically estimated on a uniform patch. It is
important to note that only the PSD of the reference image and not the
reference image itself is needed for the acutance computation. For this
reason, this method is referred to as Reduced-Reference (RR) acutance in the
rest of this paper. However, this approach is hindered by the denoising
algorithms integrated into cameras. Not only these algorithms interfere with
the noise PSD estimation but also they behave differently in uniform and
textured regions. In this context, Kirk et al. [11] propose to compute the MTF
using the cross-power spectral density of the resulting and the reference
images. This method assumes an effective registration of the chart. Sumner et
al. [12] then modified Kirk’s computation in order to make it more robust to
slight misregistration. Since this method takes fully advantage of the
reference image, it is referred to as Full-Reference (FR) acutance in the rest
of this paper.

In conclusion, state-of-the-art techniques typically allow to obtain a good
estimation of devices’ MTF and then of the acutance. However, it has been
shown that acutance itself does not always reflect very well the human quality
judgment. This observation calls for learning-based methods that aim at
reproducing the score of human experts evaluating the images.

###  II-B Learning-based methods

In opposition to MTF-based techniques described in the previous section,
learning-based methods require annotated datasets. Early datasets (LIVE [13],
CSIQ [14], TID2008 [15] and its extension TID2013 [16]) consist of noise-free
images and subjective preference scores for different types and levels of
artificially introduced distortions. These distortions mostly correspond to
compression or transmission scenarios. Only some of them, such as intensity
shift, additive Gaussian noise, contrast changes or chromatic aberrations, are
relevant to the problem of camera evaluation, and those distortions do not
encompass the complexity of real systems. Conversely, the KonIQ10k [17]
dataset consists of samples from a larger public media database with unknown
distortions claimed to be authentic. Similarly, the LIVE In the Wild [18]
database consists of 1161 images shot with various devices with unique content
for the images. For both datasets, crowd-sourced studies were conducted in
order to collect opinions about the quality of these images. Interestingly,
the size of these datasets is sufficiently large to allow the use of deep
learning approaches. In particular, Varga et al. [19] propose to fine-tune a
ResNet-101 network [20] equipped of spatial pyramid pooling layers [21] in
order to accept the full images at a non-fixed size to predict the user
preference score. This deep approach leads to the best performances on these
authentic distortions benchmarks.

Despite the good results of this method, it cannot be directly applied to our
problem. One of the reasons is that we aim at estimating the quality of the
full image produced by the camera: resizing it would alter the fine details we
try to estimate, and without resizing we undergo memory usage issues. To
tackle this problem, we propose to work on patches randomly cropped in the
input image.

More recently, Yu et al. [22] collected a dataset of 12 853 natural photos
from Flickr and used Amazon Mechanical Turk to rate them according to image
quality defects: exposition, white balance, color saturation, noise, haze,
undesired blur, composition. They used a multi-column CNN: Using GoogLeNet
[23], the first layers are shared, and then it is separated according to each
attribute. For learning efficiency purposes, they frame the problem of each
defect regression as an 11-class classification task.

Nevertheless, these deep learning approaches [21, 22] tackle a slightly
different problem than ours. They are designed to evaluate image quality
attributes for any input image while we address the problem of evaluating
devices using a known common chart.

To address the limitations of both MTF-based and learning-based methods, we
now introduce our novel deep regression framework for texture quality
estimation.

Figure 1: Pipeline overview of DR2S: Our method is divided into three main
stages. First, we perform a naive training using patches randomly cropped over
the chart. Second, we compute a map that indicates discriminant regions.
Third, we perform a final training using only a selected region.

##  III Method

###  III-A DR2S Overview

In this section, we detail the proposed method for estimating texture quality.
We formulate this task as a regression problem. We assume that we dispose of a
training dataset composed of N color images $(X_{1},\ldots,X_{N})$ of
dimensions $H\times W$ with the corresponding texture quality scores
$(Y_{1},\ldots,Y_{N})\in\mathbb{R}^{N}$. Since we are interested in estimating
the camera quality, each score $Y_{n}$ corresponds to a quality score for one
device at a specific lighting condition. Note that, several training images
can be taken with the same device. Importantly, we aim at computing quality
scores that coincide with human judgment, the ground-truth texture quality
scores are provided by human annotators (See Sec.IV-A for more details). We
aim at training a ConvNet $\phi(\cdot,\theta)$ with parameters $\theta$. We
consider that these images depict the same chart and are taken with different
cameras and different lighting conditions.

Our method is based on the observation that all image regions are not equally
suited to predict the overall image quality. For example, regions with uniform
textures will be rendered similarly by any device independently of its
quality. Conversely, other regions with rich and fine texture details are much
more discriminating since they will be differently captured by different
devices. Based on this observation, we propose the algorithm illustrated in
Fig. 1. The method is divided into three stages: First, we train a first
neural network ignoring this problem of unsuited regions. Second, we employ
this network for estimating which image regions are the most suitable for
measuring image quality. Finally, the network is re-trained using only
selected regions. In the following we provide the details of the three stages
of our algorithm.

###  III-B Initial Training

The goal of this first stage is to train a neural network that can be used to
identify relevant image regions. To this end, we propose to train a network to
regress the quality score from an input image patch. This initial network is
later used in the second stage of our pipeline to identify discriminant
regions (See Sec.III-C). We train this deep convNet on random crops extracted
from the training images $\\{X_{1},\ldots,X_{N}\\}$. These crops are randomly
selected across the images with a uniform distribution. In all our
experiments, we use the widely used Resnet-50 network pre-trained on ImageNet
[24] where the final classification layer is replaced by a linear fully
connected layer. Since some patches are not discriminant, this initial
training suffers from instability and optimization issues. Consequently, we
use the Huber loss that reports better performance in the presence of noisy
samples [25]:

$\mathcal{L}(y,\hat{y})=\left\\{\begin{array}[]{ll}\frac{1}{2}(y-\hat{y})^{2}&\text{
for }|y-\hat{y}|\leq\delta\\\ \delta|y-\hat{y}|-\frac{1}{2}\delta^{2}&\text{
otherwise }\end{array}\right.$ (1)

where $y$ and $\hat{y}$ denote the annotated and predicted scores, and
$\delta\in\mathbb{R^{+}}$ is a threshold. At the end of this stage, we obtain
a network that estimates the quality of a device from an input patch.

###  III-C Region Selection

In the second stage of our pipeline, we use our previously trained network to
predict image regions that are most discriminating for quality measurement. We
produce a map $M\in[0,1]^{H\times W}$ that indicates the relevance of each
location of the chart to estimate texture quality. This map will allow us to
select a suitable region to train our convNet in the last stage of our
pipeline.

In order to estimate a single map $M$ for all the training images, we first
register all the training images. We employ the following algorithm to align
the images on the image with the highest resolution (47 MP). First, we detect
points of interest. Then, we extract local AKAZE descriptors [26] and,
finally, we estimate a homography for every image. Image warping is
implemented using bicubic sampling. Note that, while the map computation
requires this warping alignment step that may affect the performance, training
and prediction can be performed on the original images. We now assume that the
training images $\\{X_{1},\ldots,X_{N}\\}$ are registered.

To estimate the map $M$, we propose to use the network $\phi$ trained in the
first stage. Let $\Psi\in\mathbb{R}^{H^{\prime},W^{\prime},C}$ be the feature
tensor outputted by the backbone network for a given input image. In our case,
since we employ a ResNet-50 network, $\Psi$ corresponds to the tensor before
the Global Average Pooling (GAP) layer. Since Resnet-50 is a fully
convolutional network, the dimension $H^{\prime}$ and $W^{\prime}$ depends on
the input image dimensions, while the number of channels is fixed (i.e.,
$C=2048$). Let $A\in\mathbb{R}^{C}$ and $b\in\mathbb{R}$ be the trained
parameters of the final regression layer obtained in the first stage. The
network prediction is given by:

$\hat{y}=\sigma\left(A^{\top}\cdot\text{GAP}(\Psi)+b\right)$ (2)

where $\sigma$ denotes the sigmoid function. While the network returns one
single output scalar per input image, we want to obtain one value per pixel
location. In order to adapt the class activation map framework [27] to our
regression setting, we propose to compute a score for every feature map
location $(h,w)\in[1:H^{\prime}]\times[1:W^{\prime}]$:

$S_{h,w}=\sigma\left(A^{\top}\cdot\Psi[h,w])+b\right)$ (3)

where $\Psi[h,w]\in\mathbb{R}^{C}$ denotes the feature vector at the location
$(h,w)$. Note that, the resulting map
$S=(S_{h,w})_{(h,w)\in[1:H^{\prime}]\times[1:W^{\prime}]}$ has a dimension
$H^{\prime}\times W^{\prime}$ that is different from the initial input image
size $H\times W$. This size difference depends on the network architecture. In
the case of Resnet-50, we obtain a ratio 32 between the input and the feature
map dimensions. Therefore, we resize the score map $S$ to the dimension
$H\times W$ using bicubic-sampling. This procedure, is applied to every image
of the training set. Thus, we obtain the set of score maps
$\\{\tilde{S}_{1},\ldots,\tilde{S}_{N}\\}$

We propose to define the confidence score map $M$ as the location-wise
variance of the score maps $\tilde{S}_{n}$ over the whole training set. The
motivation for this choice is that, discriminating regions have a higher
variance than non-discriminating ones. Indeed, we observe that the scores
produced by the ConvNet $\phi(\cdot,\theta)$ over non-discriminative regions
tend to have small variance. Conversely, on discriminating patches, the
networks predict values with a wide range leading to high variance.

###  III-D Final Training and Prediction

In the last stage of our pipeline, we select the chart region with the highest
confidence score value in $M$. In our preliminary experiments, we observed
that using a region width approximately six times larger than the network
input size leads to good performance. In our case, we use a square region of
$1200\times 1200$ pixels. In this region, we select random patches that are
used as a training set. We re-train the network $\phi$, starting again from
ImageNet pre-training weights.

At test time, assuming an image with an unknown quality score, we extract
patches in the selected region. The final score is given by the average over
the different patches.

##  IV Experiments

In this section, we perform a thorough experimental evaluation of the proposed
pipeline. We implemented our method using Tensorflow and Keras. When training
the ConvNets (stages 1 and 3), we employ Adam optimizer following [28], with a
starting learning rate of $3\cdot 10^{-3}$ with a decay of 0.1 every 40 epochs
for a total of 120 epochs. In our model, we assume that all the images have
the same resolution. The reason for this choice is that we want the image
details to be analyzed at the same scale, as a human observer would do. In
practice, the resolution depends on the device. Therefore, we preprocess all
the training images resizing them to the highest resolution of the dataset
using bicubic upsampling. This solution is preferred to downsampling to a
common lower resolution since texture quality is not invariant to
downsampling. In addition, due to possible lens shading, we remove the sides
of the images.

###  IV-A Datasets

##### Charts and devices

As there is no well-established reference dataset for our problem, we
collected annotated data using three different charts.

  * •

Still-Life: First, we use the chart displayed in Fig. 3(a). This dataset is
referred to as _Still-Life_. The chart is designed to evaluate several image
quality attributes and to present diverse content: Vivid colors for color
rendering, fine details, uniform zones, portraits as well as resolution lines
and a low-quality Dead-Leaves version. Images are acquired using 140 different
smartphones and cameras from different brands commonly available in the
consumer market. In Fig. 2, we provide an example patch captured using three
different cameras. The left image corresponds to a high-quality device while
the two others are obtained with low-quality devices. It illustrates the
nature of distortions that appear in this dataset with different intensity. To
obtain a larger database and predictions robust to lighting conditions, we
shoot the chart using five different lighting conditions: 5 lux tungsten, 20
lux tungsten, 100 lux tungsten, 300 lux TL84, 1000 lux D65. Note that process
is repeated for every device.

  * •

Gray-DL: Second, we employ the dead-leave chart proposed in [10]. As mentioned
in Sec II-A, this chart depicts gray-scale circles with random radius and
locations. In all our experiments, we refer to this dataset as _Gray-DL_. We
use the same five lighting conditions and devices as for the _Still-Life_
chart.

  * •

Color-DL: Finally, we complete our experiment using the dead-leaves chart
proposed in [12]. By opposition to _Gray-DL_ this chart is colored (an image
can be found at [29]). For this chart, we employed a limited number of
devices. More precisely, we employ only 14 devices with the same five lighting
conditions as for the other charts. The low number of devices is especially
challenging for learning-based approaches.

Figure 2: Patches of high (left) and low-quality (center and right) images
from our _Still-Life_ dataset

(a) _Still-Life_ Chart

(b) _Gray-DL_ Chart

Figure 3: Still-Life Chart used in our experiments. The _Still-Life_ chart
contains many diverse objects with varying colors and textures while the
_Gray-DL_ chart depicts random gray-scale circles.

##### Annotations

In order to train and evaluate the different methods we need to provide
ground-truth annotation for each device. Note that, the annotation must be
provided for each pair of device and lighting conditions. To obtain quality
annotations that are a reliable proxy of the perceived visual quality,
annotations are provided by human experts. Images to be evaluated were
inserted among a fixed set of 42 references, and very high-quality prints were
provided to help them judge of the authenticity of the details. Annotators
were asked to compare the images using the same field of view for every image,
using calibrated high-quality monitors where images are displayed without
applying any down-sampling but with a possible digital zoom for the lower
resolution image. Each position among the set of references is assigned a
score between 0 and 1. In the case of the Dead-Leaves charts, since the charts
are unnatural images, human perceptual annotation is problematic. Therefore,
we chose to use the annotations obtained on the _Still-Life_ also for the
dead-leaves charts, rather than annotating the images. The _Still-Life_ chart
contains diverse textures similar to what real images would contain. In this
way, we obtain a subjective device evaluation in a setting more similar to
real-life scenarios.

###  IV-B Metrics

In our problem, relying on standard classification or regression metrics is
not straightforward. Indeed, MTF-based methods predict a quality score that is
not directly comparable to the score provided by human annotators. A
straightforward alternative could consist in computing the correlation between
the predictions and the annotation. However, the underlying assumption that
the predictions of each method correlate linearly with our annotations may not
hold and bias the evaluation. Therefore, we decided to rely on two distinct
metrics based on the correlation of the rank-order. First, we adopt the
Spearman Rank-Order Correlation Coefficient (_SROCC_) defined as the linear
correlation coefficient of the ranks of predictions and annotations. Second,
we report the Kendall Rank-Order Correlation Coefficient (_KROCC_) defined by
the difference between concordant and discordant pairs divided by the number
of possible pairs. The key advantage of this second metric lies in its
robustness to outliers.

For all visual charts, the dataset is split into training and test sets as
follows. First, among the devices we use in our experiments several are
produced by the same brand. So, to avoid bias between training and test, we
impose no brand-overlap between training and test sets. Second, as a
consequence of such constrain, a limited number of brands may appear in the
test set. To avoid evaluation biases towards specific brands, we use a
$k$-fold cross-validation with $k=16$.

In order to measure the impact of the number of devices on the performance, we
perform experiments with a variable number of devices. For all the experiments
on the _Gray-DL_ and _Still-life_ charts, we report the results obtained using
subsets of size 20, 60, 100 and 140 devices. For a given number of devices,
each experience is performed over the same devices set. Note that, for every
method, the complete pipeline is repeated independently for every subset.

###  IV-C Ablation study

In order to experimentally justify our proposed method, we compare three
different versions of our model:

  * •

_Random Patch_ : In this approach, random patches are selected from the whole
chart at both training and testing time.

  * •

_Random Region_ : We then restrict the random patch extraction to a single
zone, chosen randomly. We report the average over five random regions.

  * •

_Selected Region_ : In this model we employ our full pipeline as described in
Sec. III. In particular, training and test are performed using the selected
region.

In these three models, we employ a ResNet-50 backbone trained using the same
optimization hyper-parameters.

TABLE I: Ablation Study: we measure the impact of region selection comparing three baseline models on the _Still-life_ chart. _SROCC_ and _KROCC_ metrics are reported. Number of devices | 20 | 60 | 100 | 140 |  | 20 | 60 | 100 | 140
---|---|---|---|---|---|---|---|---|---
| SROCC |  | KROCC
Random Patch | 0.626 | 0.818 | 0.784 | 0.806 |  | 0.433 | 0.617 | 0.588 | 0.613
Random Region | 0.795 | 0.863 | 0.866 | 0.879 |  | 0.606 | 0.680 | 0.682 | 0.700
Selected Region (Full model) | 0.830 | 0.912 | 0.890 | 0.900 |  | 0.638 | 0.740 | 0.716 | 0.728

The results obtained on the _Still-life_ chart are reported in Table.I. First,
when using the _random patches_ variant, the model trained on 20 devices
performs poorly both in terms of _SROCC_ and _KROCC_ compared to other
variants. In this case, we see that it is required to dispose of at least 60
devices to get satisfying performances. Second, we observe that restricting
random patches extraction to a region randomly selected leads to better
performance than if we do not restrict to this region. The gain is visible for
every number of devices and for both metrics. It may be explained that the
decreased diversity in content leads to a ConvNet that is specialized in a
specific region of the chart. In other words, the benefit of a more restrained
input diversity is larger than the benefit of a larger and more diverse
training set. Finally, our full model reaches the best performance for both
metrics and for every number of devices. This better performance independently
of the training sub-set demonstrates the robustness of the proposed method.
Interestingly, we obtain performances with 20 devices similar to the
performance of the _Random Patch_ model with 140 devices. Overall, this
ablation study illustrates the benefit of selecting specific regions for
texture quality measurement.

###  IV-D Qualitative analysis of our region selection

In order to further study the outcome of our region selection algorithm, we
display the resulting map (Fig. 4) of relevant zones.

Figure 4: Normalized discriminant-region map $S$ (better viewed with a digital
zoom). For display, we employ histogram equalization for normalization and
obtain values from 0 to 1.

We observe that uniform regions are considered by our algorithm as the least
discriminant for texture quality assessment. In particular, this is visible in
the bottom-right regions on the black square patch. On the contrary, regions
with low contrast and many small details appear to be more discriminant (see
around the banknote region). Results on wooden regions seem to depend on wood
grain.

This analysis is performed considering all the images. We now propose to
analyze the regions that discriminate devices among only low quality or only
high-quality images. For this analysis, the test set is split to according to
the ground-truth score. In this way, we compute two discriminant maps. Two
small crops of these two maps are shown in Fig.5. Interestingly, we observe
restricting our analysis to high quality or lower quality images leads to
differences in results. For example, we observe that the resolution lines (in
the bottom row of 5) discriminate for low-quality images, but not for higher
quality images. Conversely, areas exhibiting only very fine details are not
the most useful for low-quality images. In particular, the forehead of the man
is not discriminant among low-quality images, while this region is highly
discriminant among high-quality images. It shows that the region with very
fine details are discriminant only among high-quality images since these
details are completely distorted by all the low-quality devices.

(a) Low quality images

(b) High quality images

Figure 5: Comparison of discriminant-region maps for high-quality images and
low quality. We display two patches extracted from the confidence maps
obtained when using texture quality maps only from high (ie. Left) and low
(ie. Right) quality images.

###  IV-E Comparison to state of the art

TABLE II: Comparison of deep learning systems on different charts to [10]._SROCC_ and _KROCC_ metrics are reported. Number of devices |  | 20 | 60 | 100 | 140 |  | 20 | 60 | 100 | 140
---|---|---|---|---|---|---|---|---|---|---
Method | Chart |  | SROCC |  | KROCC
RR Acutance [10] | _Gray-DL_ |  | 0.704 | 0.794 | 0.747 | 0.788 |  | 0.533 | 0.595 | 0.592 | 0.592
ResNet [20] | _Gray-DL_ |  | 0.641 | 0.795 | 0.792 | 0.824 |  | 0.464 | 0.598 | 0.592 | 0.630
$DR^{2}S$ (Ours) | _Still-Life_ |  | 0.830 | 0.912 | 0.890 | 0.900 |  | 0.638 | 0.740 | 0.716 | 0.728

In this section, we compare the performance of our approach to existing
methods. This comparison is twofold since both the methods and the charts need
to be compared. We perform two sets of experiments. In the first set of
experiments, we compare different methods on the two large datasets recorded
with the _Gray-DL_ and the _Still-Life_ charts and the same 140 devices. The
second set of experiments consists of a comparison of the devices on the
_Color-DL_ chart. This second set of experiments is highly challenging for
learning-based methods because of the limited amount of training data.

#### Large database experiments

First, in our preliminary experiments, we observed that, for this experience,
adding a small amount of Gaussian noise and random change in exposition leads
to better performance. This data-augmentation is performed on the fly on every
training patch. Our main competitor is the _RR acutance_ methods proposed in
[10]. For the acutance computation, viewing conditions were set to 120
centimeters printing height and 100 centimeters viewing distance. Note that,
the _RR acutance_ method is intrinsically designed for the Dead-Leaves charts
and cannot be used for the _Still-Life_ chart. We include a second deep
learning-based method for the _Gray-DL_ chart in our comparison. This approach
consists of a ResNet-50 [20] where the classification layer is replaced by 3
additional fully-connected layers and a linear regression layer. For this
approach, we employ the _Random patch_ strategy described in Sec. IV-C inside
of the texture region. Importantly, we do not report the performance of _DR
2S_ on the _Gray-DL_ chart since the chart is designed to be uniformly
discriminant for texture quality assessment.

Quantitative results are reported in Table. II. First, we observe that with a
limited number of devices for training (e.g. 20 devices), _RR Acutance_
performs better than ResNet-50. However, the proposed approach clearly
outperforms the texture-MTF based method ($+0.126$ and $+0.105$ in _SROCC_ and
_KROCC_ , respectively). It shows that when few training samples are
available, selecting the appropriate regions is essential for good
performance. ResNet performance increases with the number of devices: with 140
devices, ResNet-50 clearly outperforms _RR Acutance_ according to both metrics
showing the potential of learning-based methods. While comparisons between
results obtained using different charts must be interpreted with care, this
result clearly shows that a learning-based approach can be intrinsically
better than acutance-based methods using the exact same input images. Finally
our _DR 2S_ method on the _Still-life_ chart leads to the best results
according to both metrics and for every number of devices.

#### Small database experiments

Concerning the second set of experiments, we compared the different methods on
the _Color-DL_ chart. Note that the 14 devices of the _Color-DL_ chart are a
subset of the devices of the _Gray-DL_ and _Still-Life_ charts. Consequently,
we also performed experiments using only these 14 devices on these two other
charts. Note that this setting is very challenging for the two learning-based
methods (ResNet and _DR 2S_) because of the limited amount of training data.
Therefore, for the two learning methods, we perform two experiments. First,
training and test are performed using 14-fold cross-validation on the exact
same data as the other methods. Second, we train our model on the complete
database and test on the 14 devices in common with the _Color-DL_ chart. These
two variants are referred to as _Restricted_ and _Full_. Again, we do not
report the performance of _DR 2S_ on the _Gray-DL_ and _Color-DL_ charts since
the chart is designed to be uniformly discriminant. Results are reported in
Table III.

TABLE III: State of the art comparison : Performance on the 14 devices database. Deep learning systems on perceptual and Gray-DL are compared to [10] and [12] Method | Chart |  | SROCC | KROCC
---|---|---|---|---
_FR Acutance_ [12] | _Color-DL_ |  | 0.701 | 0.544
_RR Acutance_ [10] | _Gray-DL_ |  | 0.714 | 0.552
ResNet - Restricted | _Gray-DL_ |  | 0.640 | 0.463
ResNet - Full | _Gray-DL_ |  | 0.780 | 0.598
$DR^{2}S$ \- Restricted | _Still-Life_ |  | 0.746 | 0.569
$DR^{2}S$ \- Full | _Still-Life_ |  | 0.873 | 0.702

First, we observe that the MTF-based methods perform similarly on the color
and gray-scale dead leave charts. It shows that better performance of the
proposed model on the _Still-Life_ chart is not due to the lack of colors in
_Gray-DL_ but to its content. Second, using the restricted database, both
learning-based methods, ResNet and $DR^{2}S$, under-perform MTF-based
predictions. However, when the amount of training data is sufficient, both
methods outperform _FR Acutance_ and _RR Acutance_.

##  V Conclusion and Future Work

In this paper, we proposed DR2S, a method which learns to estimate a
perceptual quality score. To this end, our algorithm selects the chart region
that is the most suitable for texture quality assessment. Our results also
suggest that, if enough training samples are available, learning-based methods
outperform MTF-based methods. A limitation of our method is that we select
only a single region. However, texture quality is known to be multi-
dimensional. Consequently, as future work, we plan to extend our method to
multiple regions in order to highlight several complementary discriminant
features and better measure the intrinsic qualities of a device.

## References

  * [1] M. C. Trinidad, R. M. Brualla, F. Kainz, and J. Kontkanen, “Multi-view image fusion,” in _Proceedings of the IEEE International Conference on Computer Vision_ , 2019.
  * [2] S. W. Hasinoff, D. Sharlet, R. Geiss, A. Adams, J. T. Barron, F. Kainz, J. Chen, and M. Levoy, “Burst photography for high dynamic range and low-light imaging on mobile cameras,” _ACM Transactions on Graphics (TOG)_ , 2016.
  * [3] M. Čadík, M. Wimmer, L. Neumann, and A. Artusi, “Image attributes and quality for evaluation of tone mapping operators,” in _National Taiwan University_. Citeseer, 2006.
  * [4] O. van Zwanenberg, S. Triantaphillidou, R. Jenkin, and A. Psarrou, “Edge detection techniques for quantifying spatial imaging system performance and image quality,” in _Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition Workshops_ , 2019.
  * [5] R. Zhang, P. Isola, A. A. Efros, E. Shechtman, and O. Wang, “The unreasonable effectiveness of deep features as a perceptual metric,” in _Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition_ , 2018.
  * [6] G. Boreman, _Modulation Transfer Function in Optical and Electro-optical Systems_. Society of Photo Optical, 2001\.
  * [7] C. Loebich, D. Wueller, B. Klingen, and A. Jaeger, “Digital camera resolution measurements using sinusoidal siemens stars,” in _Digital Photography III_. International Society for Optics and Photonics, 2007.
  * [8] G. Matheron, _Random sets and integral geometry_ , ser. Wiley series in probability and mathematical statistics: Probability and mathematical statistics. Wiley, 1975.
  * [9] Y. Gousseau and F. Roueff, “Modeling occlusion and scaling in natural images,” _Multiscale Modeling & Simulation_, 2007.
  * [10] F. Cao, F. Guichard, and H. Hornung, “Measuring texture sharpness of a digital camera,” in _Digital Photography V_. International Society for Optics and Photonics, 2009.
  * [11] L. Kirk, P. Herzer, U. Artmann, and D. Kunz, “Description of texture loss using the dead leaves target: current issues and a new intrinsic approach,” in _Digital Photography X_. International Society for Optics and Photonics, 2014.
  * [12] R. C. Sumner, R. Burada, and N. Kram, “The effects of misregistration on the dead leaves crosscorrelation texture blur analysis,” _Electronic Imaging_ , 2017.
  * [13] H. R. Sheikh, M. F. Sabir, and A. C. Bovik, “A statistical evaluation of recent full reference image quality assessment algorithms,” _IEEE Transactions on image processing_ , 2006.
  * [14] E. C. Larson and D. M. Chandler, “Most apparent distortion: full-reference image quality assessment and the role of strategy,” _Journal of electronic imaging_ , 2010.
  * [15] N. Ponomarenko, V. Lukin, A. Zelensky, K. Egiazarian, M. Carli, and F. Battisti, “Tid2008-a database for evaluation of full-reference visual quality assessment metrics,” _Advances of Modern Radioelectronics_ , 2009\.
  * [16] N. Ponomarenko, L. Jin, O. Ieremeiev, V. Lukin, K. Egiazarian, J. Astola, B. Vozel, K. Chehdi, M. Carli, F. Battisti _et al._ , “Image database tid2013: Peculiarities, results and perspectives,” _Signal Processing: Image Communication_ , 2015.
  * [17] V. Hosu, H. Lin, T. Sziranyi, and D. Saupe, “Koniq-10k: An ecologically valid database for deep learning of blind image quality assessment,” _IEEE Transactions on Image Processing_ , 2020.
  * [18] D. Ghadiyaram and A. C. Bovik, “Massive online crowdsourced study of subjective and objective picture quality,” _IEEE Transactions on Image Processing_ , vol. 25, no. 1, pp. 372–387, 2016.
  * [19] D. Varga, D. Saupe, and T. Szirányi, “Deeprn: A content preserving deep architecture for blind image quality assessment,” in _2018 IEEE International Conference on Multimedia and Expo (ICME)_. IEEE, 2018.
  * [20] K. He, X. Zhang, S. Ren, and J. Sun, “Deep residual learning for image recognition,” in _Proceedings of the IEEE conference on computer vision and pattern recognition_ , 2016.
  * [21] ——, “Spatial pyramid pooling in deep convolutional networks for visual recognition,” _IEEE transactions on pattern analysis and machine intelligence_ , 2015.
  * [22] N. Yu, X. Shen, Z. Lin, R. Mech, and C. Barnes, “Learning to detect multiple photographic defects,” in _2018 IEEE Winter Conference on Applications of Computer Vision (WACV)_. IEEE, 2018\.
  * [23] C. Szegedy, W. Liu, Y. Jia, P. Sermanet, S. Reed, D. Anguelov, D. Erhan, V. Vanhoucke, and A. Rabinovich, “Going deeper with convolutions,” in _Proceedings of the IEEE conference on computer vision and pattern recognition_ , 2015.
  * [24] J. Deng, W. Dong, R. Socher, L.-J. Li, K. Li, and L. Fei-Fei, “Imagenet: A large-scale hierarchical image database,” in _2009 IEEE conference on Computer Vision and Pattern Recognition_ , 2009.
  * [25] M. Carvalho, B. Le Saux, P. Trouvé-Peloux, A. Almansa, and F. Champagnat, “On regression losses for deep depth estimation,” in _2018 25th IEEE International Conference on Image Processing (ICIP)_. IEEE, 2018.
  * [26] P. F. Alcantarilla and T. Solutions, “Fast explicit diffusion for accelerated features in nonlinear scale spaces,” _IEEE Trans. Pattern. Anal. Mach. Intell_ , 2011.
  * [27] B. Zhou, A. Khosla, A. Lapedriza, A. Oliva, and A. Torralba, “Learning deep features for discriminative localization,” in _Proceedings of the IEEE conference on computer vision and pattern recognition_ , 2016.
  * [28] S. Lathuilière, P. Mesejo, X. Alameda-Pineda, and R. Horaud, “A comprehensive analysis of deep regression,” _IEEE TPAMI_ , 2020.
  * [29] “Imatest website.” [Online]. Available: https://www.imatest.com/docs/random-cross


a_{i_{n}}b_{i_{n}}^{\prime})=\operatorname{tp}(a_{j_{1}}b_{j_{1}}^{\prime}\ldots
a_{j_{n}}b_{j_{n}}^{\prime}).$

So $\operatorname{tp}(a_{i_{1}}[b_{i_{1}}]\ldots
a_{i_{n}}[b_{i_{n}}])\subseteq\operatorname{tp}(a_{j_{1}}[b_{j_{1}}]\ldots
a_{j_{n}}[b_{j_{n}}])$, and the claim follows by maximality of types. ∎

######  Theorem 10.17.

thm:hyperimaginaries-preserving-hausdorff-thickness The following properties
of $T$ are preserved when adding hyperimaginaries:

  * •

Hausdorff,

  * •

semi-Hausdorff,

  * •

thick.

That is, if $T$ has the property then $T^{\mathcal{E}}$ has it as well.

###### Proof.

Hausdorff. Let $a[b]\not\equiv a^{\prime}[b^{\prime}]$. Then there is
$\varphi\in\operatorname{tp}(a[b])$ such that
$\varphi\not\in\operatorname{tp}(a^{\prime}[b^{\prime}])$. So there is a
negation $\psi\in\operatorname{tp}(a^{\prime}[b^{\prime}])$ of $\varphi$. By
lem:translation-lemma we have $\Sigma_{\varphi}$ and $\Sigma_{\psi}$ are
consistent while $\Sigma_{\varphi}\cup\Sigma_{\psi}$ is inconsistent.

Fix some type $q$ of $T$ such that $\Sigma_{\psi}\subseteq q$. We will produce
formulas $\alpha_{q}$ and $\beta_{q}$ such that
$\Sigma_{\varphi}\cup\\{\alpha_{q}\\}$ is inconsistent, $\beta_{q}\not\in q$
and $T\models\forall xy_{r}(\alpha_{q}(x,y_{r})\vee\beta_{q}(x,y_{r}))$. Let
$p\supseteq\Sigma_{\varphi}$ be a type of $T$. Then because $T$ is Hausdorff
there are formulas $\chi_{p}$ and $\theta_{p}$ such that $\chi_{p}\not\in p$
and $\theta_{p}\not\in q$, while $T\models\forall
xy_{r}(\chi_{p}(x,y_{r})\vee\theta_{p}(x,y_{r}))$. Then
$\Sigma_{\varphi}\cup\\{\chi_{p}:p\supseteq\Sigma_{\varphi}\\}$ is
inconsistent, so there are $p_{1},\ldots,p_{n}$ such that
$\Sigma_{\varphi}\cup\\{\chi_{p_{1}}\wedge\ldots\wedge\chi_{p_{n}}\\}$ is
inconsistent. We can now take $\alpha_{q}$ to be
$\chi_{p_{1}}\wedge\ldots\wedge\chi_{p_{n}}$ and $\beta_{q}$ to be
$\theta_{p_{1}}\vee\ldots\vee\theta_{p_{n}}$.

Now $\Sigma_{\psi}\cup\\{\beta_{q}:q\supseteq\Sigma_{\psi}\\}$ is
inconsistent. So there are $q_{1},\ldots,q_{k}$ such that
$\Sigma_{\psi}\cup\\{\beta_{q_{1}}\wedge\ldots\wedge\beta_{q_{k}}\\}$ is
inconsistent. We set $\beta=\beta_{q_{1}}\wedge\ldots\wedge\beta_{q_{k}}$ and
$\alpha=\alpha_{q_{1}}\vee\ldots\vee\alpha_{q_{n}}$. We then also have that
$\Sigma_{\varphi}\cup\\{\alpha\\}$ is inconsistent and $T\models\forall
xy_{r}(\alpha(x,y_{r})\vee\beta(x,y_{r}))$.

Now consider the formulas $R_{\alpha}(x,y)$ and $R_{\beta}(x,y)$. By
construction we have $T^{\mathcal{E}}\models\forall xy(R_{\alpha}(x,y)\vee
R_{\beta}(x,y))$. We claim that $R_{\alpha}\not\in\operatorname{tp}(a[b])$.
Suppose for a contradiction that $\mathfrak{M}^{\mathcal{E}}\models
R_{\alpha}(a,[b])$. Then there is $b^{*}$ with $[b^{*}]=[b]$ such that
$\mathfrak{M}\models\alpha(a,b^{*})$. Since
$\varphi\in\operatorname{tp}(a[b])=\operatorname{tp}(a[b^{*}])$ we also have
$\mathfrak{M}\models\Sigma_{\varphi}(a,b^{*})$, contradicting that
$\Sigma_{\varphi}\cup\\{\alpha\\}$ is inconsistent. So indeed
$R_{\alpha}\not\in\operatorname{tp}(a[b])$. Analogously we get that
$R_{\beta}\not\in\operatorname{tp}(a^{\prime}[b^{\prime}])$, which concludes
the proof that $T^{\mathcal{E}}$ is Hausdorff.

Semi-Hausdorff. Suppose that equality of $\mathcal{L}$-types is type-definable
by a partial $\mathcal{L}$-type $\Omega$. Then for a tuple $x$ of real
variables and a tuple $y$ of hyperimaginary variables, we consider the partial
$\mathcal{L}_{\mathcal{E}}$-type
$\Omega^{\mathcal{E}}(xy,x^{\prime}y^{\prime})$ that expresses the following:

$\exists
y_{r}y_{r}^{\prime}(\Xi(y_{r},y)\wedge\Xi(y_{r}^{\prime},y^{\prime})\wedge\Omega(xy_{r},x^{\prime}y_{r}^{\prime})).$

We claim that $\Omega^{\mathcal{E}}$ expresses equality of
$\mathcal{L}_{\mathcal{E}}$-types.

If
$\mathfrak{M}^{\mathcal{E}}\models\Omega^{\mathcal{E}}(a[b],a^{\prime}[b^{\prime}])$
then we find $c,c^{\prime}$ such that
$\mathfrak{M}^{\mathcal{E}}\models\Xi(c,[b])\wedge\Xi(c^{\prime},[b^{\prime}])\wedge\Omega(ac,a^{\prime}c^{\prime})$.
By lem:hyperimaginary-projection-type we have that $[c]=[b]$ and
$[c^{\prime}]=[b^{\prime}]$. Hence
$\varphi\in\operatorname{tp}(a[b])=\operatorname{tp}(a[c])$ iff
$\Sigma_{\varphi}\subseteq\operatorname{tp}(ac)=\operatorname{tp}(a^{\prime}c^{\prime})$
iff
$\varphi\in\operatorname{tp}(a^{\prime}[c^{\prime}])=\operatorname{tp}(a^{\prime}[b^{\prime}])$.
So $\operatorname{tp}(a[b])=\operatorname{tp}(a^{\prime}[b^{\prime}])$, as
required.

Conversely, if
$\operatorname{tp}(a[b])=\operatorname{tp}(a^{\prime}[b^{\prime}])$ then by
lem:type-hyperimaginaries-determined-by-type-reals we find $b^{\prime\prime}$
such that $[b^{\prime\prime}]=[b^{\prime}]$ and
$\operatorname{tp}(ab)=\operatorname{tp}(a^{\prime}b^{\prime\prime})$. Hence
$\models\Xi(b,[b])\wedge\Xi(b^{\prime\prime},[b^{\prime}])\wedge\Omega(ab,a^{\prime}b^{\prime\prime})$.

Thick. Let $\Theta$ express indiscernibility of a sequence of real tuples,
then

$\exists(y_{i,r})_{i<\omega}\left(\Theta((x_{i}y_{i,r})_{i<\omega})\wedge\bigwedge_{i<\omega}\Xi(y_{i,r},y_{i})\right)$

expresses indiscernibility of $(x_{i}y_{i})_{i<\omega}$ in $T^{\mathcal{E}}$.
Here we use that a sequence in $\mathfrak{M}^{\mathcal{E}}$ is indiscernible
if and only if there is an indiscernible sequence of real representatives, see
lem:indiscernible-sequences-real-representation. ∎

######  Theorem 10.18.

thm:hyperimaginaries-nsop1 The theory $T$ is NSOP1 if and only if
$T^{\mathcal{E}}$ is NSOP1.

The technique in the proof of thm:hyperimaginaries-nsop1 can also be applied
to other combinatorial properties, such as the order property, TP, TP2, IP,
etc. Of course, to do this, one first needs to write down a proper definition
of these properties for positive logic, such as def:sop1 for SOP1 or [HK21,
Definition 6.1] for TP2.

###### Proof.

One direction is trivial: if $T$ has a formula with SOP1, then so has
$T^{\mathcal{E}}$.

We prove the other direction: suppose that $T^{\mathcal{E}}$ has a formula
with SOP1, we will show that $T$ already has a formula with SOP1. So let
$\varphi(x,y;w,z)$ be an $\mathcal{L}_{\mathcal{E}}$-formula with SOP1. Here
$x$ and $w$ are tuples of real variables, and $y$ and $z$ are tuples of
hyperimaginary variables. Let $(a_{\eta}[b_{\eta}]:\eta\in 2^{<\omega})$ and
$\psi(w_{1},z_{1};w_{2},z_{2})$ be witnesses of SOP1. Let
$\Sigma_{\varphi}(x,y_{r};w,z_{r})$ and
$\Sigma_{\psi}(w_{1},z_{1,r};w_{2},z_{2,r})$ be as in lem:translation-lemma-
formula. Then

$\Sigma_{\psi}(w_{1},z_{1,r};w_{2},z_{2,r})\cup\Sigma_{\varphi}(x,y_{r},w_{1};z_{1,r})\cup\Sigma_{\varphi}(x,y_{r},w_{2};z_{2,r})$

is inconsistent. Hence there are finite $\varphi^{\prime}\in\Sigma_{\varphi}$
and $\psi^{\prime}\in\Sigma_{\psi}$ that are inconsistent with each other.
That is

$\displaystyle T\models\neg\exists
xy_{r}w_{1}z_{1,r}w_{2}z_{2,r}(\psi^{\prime}(w_{1},z_{1,r},w_{2},z_{2,r})\wedge\varphi^{\prime}(x,y_{r},w_{1},z_{1,r})\wedge\varphi^{\prime}(x,y_{r},w_{2},z_{2,r})).$
(4)

As usual, any variables not actually appearing in the formulas should be
ignored in the existential quantifier. We claim that $\varphi^{\prime}$ has
SOP1, which is witnessed by $(a_{\eta}b_{\eta}:\eta\in 2^{<\omega})$ and
$\psi^{\prime}$. We check the items in def:sop1.

  1. (i)

Let $\sigma\in 2^{\omega}$, then
$\\{\varphi(x,y,a_{\sigma|_{n}},[b_{\sigma|_{n}}]):n<\omega\\}$ is consistent.
So there are $c$ and $[d]$ such that
$\models\varphi(c,[d],a_{\sigma|_{n}},[b_{\sigma|_{n}}])$ for all $n<\omega$.
That is, we have $\Sigma_{\varphi}(c,d,a_{\sigma|_{n}},b_{\sigma|_{n}})$ for
all $n<\omega$. In particular
$\\{\varphi^{\prime}(x,y_{r},a_{\sigma|_{n}},b_{\sigma|_{n}}):n<\omega\\}$ is
consistent.

  2. (ii)

By construction, see (4).

  3. (iii)

Let $\eta,\nu\in\omega^{<\omega}$ such that $\eta^{\frown}0\preceq\nu$. Then
$\models\psi(a_{\eta^{\frown}1},[b_{\eta^{\frown}1}],a_{\nu},[b_{\nu}])$, so
$\models\Sigma_{\psi}(a_{\eta^{\frown}1},b_{\eta^{\frown}1},a_{\nu},b_{\nu})$
and in particular
$\models\psi^{\prime}(a_{\eta^{\frown}1},b_{\eta^{\frown}1},a_{\nu},b_{\nu})$.

∎

######  Definition 10.19.

def:existence-axiom We say that a theory satisfies the _existence axiom for
forking_ if $\operatorname{tp}(a/B)$ does not fork over $B$ for any $a$ and
$B$.

######  Theorem 10.20.

thm:hyperimaginaries-existence-axiom-for-forking The theory $T$ satisfies the
existence axiom for forking if and only if $T^{\mathcal{E}}$ satisfies the
existence axiom for forking.

###### Proof.

One direction is immediate: anything witnessing forking in $T$ will also be in
$T^{\mathcal{E}}$. We prove the other direction. So assume that there is
$\operatorname{tp}(a[b]/C[D])$ that forks over $C[D]$. That is, it implies a
(possibly infinite) disjunction $\bigvee_{i\in I}\varphi_{i}(xy,e^{i}[f^{i}])$
with $\varphi_{i}(xy,e^{i}[f^{i}])$ dividing over $C[D]$ for each $i\in I$.
For each $i\in I$ we let $(e^{i}_{j}[f^{i}_{j}])_{j\in J}$ be a long enough
$C[D]$-indiscernible sequence with $e^{i}_{0}[f^{i}_{0}]=e^{i}[f^{i}]$ such
that $\\{\varphi_{i}(xy,e^{i}_{j}[f^{i}_{j}]):j\in J\\}$ is inconsistent. By
lem:type-hyperimaginaries-determined-by-type-reals we may assume that
$e^{i}_{j}f^{i}_{j}\equiv e^{i}f^{i}$ for every $j\in J$. We claim that
$\Sigma_{\varphi_{i}}(x,y_{r},e^{i},f^{i})$ (see lem:translation-lemma-
formula) divides over $CD$ for all $i\in I$. Note that $\Sigma_{\varphi_{i}}$
may contain parameters from $CD$.

To prove the claim let $k$ be such that
$\\{\varphi_{i}(xy,e^{i}_{j}[f^{i}_{j}]):j\in J_{0}\\}$ is inconsistent for
all $J_{0}\subseteq J$ with $|J_{0}|=k$. So $\bigcup_{j\in
J_{0}}\Sigma_{\varphi_{i}}(x,y_{r},e^{i}_{j},f^{i}_{j})$ is inconsistent for
all such $J_{0}$. Let $(e_{n}f_{n})_{n<\omega}$ be a $CD$-indiscernible
sequence based on $(e^{i}_{j}f^{i}_{j})_{j\in J}$ over $CD$. Then there are
$j_{1}<\ldots<j_{k}\in J$ such that $e_{1}f_{1}\ldots
e_{k}f_{k}\equiv_{CD}e^{i}_{j_{1}}f^{i}_{j_{1}}\ldots
e^{i}_{j_{k}}f^{i}_{j_{k}}$, so
$\bigcup_{n<\omega}\Sigma_{\varphi_{i}}(x,y_{r},e_{n},f_{n})$ is inconsistent.
We conclude that $\Sigma_{\varphi_{i}}(x,y_{r},e^{i},f^{i})$ divides over
$CD$, as claimed.

By the claim there is $\psi_{i}(x,y_{r},e^{i},f^{i})$ that is implied by
$\Sigma_{\varphi_{i}}(x,y_{r},e^{i},f^{i})$ such that
$\psi_{i}(x,y_{r},e^{i},f^{i})$ divides over $CD$, for all $i\in I$. Let
$p=\operatorname{tp}(a[b]C[D])$, then $\Sigma_{p}(x,y_{r},C,D)$ implies
$\bigvee_{i\in I}\Sigma_{\varphi_{i}}(x,y_{r},e^{i},f^{i})$. We thus have that
$\Sigma_{p}(x,y_{r},C,D)$ implies $\bigvee_{i\in
I}\psi_{i}(x,y_{r},e^{i},f^{i})$. So $\Sigma_{p}(x,y_{r},C,D)$ forks over
$CD$. ∎

In the discussion following Definition 4.1 in [Kim] it is stated that one may
produce results for Kim-independence for the hyperimaginary extension
$M^{\text{heq}}$ of a first-order structure $M$ parallel with those for first-
order structures, provided that $M^{\text{heq}}$ satisfies the existence axiom
for forking (which, by the above theorem, is equivalent to the assumption that
$T$ satisfies this axiom). More generally, one can ask if our results on Kim-
independence over models in thick NSOP1 theories can be extended to arbitrary
base sets assuming the existence axiom for forking:

######  Question 10.21.

q:sets Suppose $T$ is a thick positive NSOP1 theory satisfying the existence
axiom for forking. Can $\mathop{\mathchoice{\displaystyle\kern 5.71527pt\hbox
to0.0pt{\hss$\displaystyle\mid$\hss}\lower 3.87495pt\hbox
to0.0pt{\hss$\displaystyle\smile$\hss}\kern 5.71527pt}{\textstyle\kern
5.71527pt\hbox to0.0pt{\hss$\textstyle\mid$\hss}\lower 3.87495pt\hbox
to0.0pt{\hss$\textstyle\smile$\hss}\kern 5.71527pt}{\scriptstyle\kern
2.80048pt\hbox to0.0pt{\hss$\scriptstyle\mid$\hss}\lower 1.89871pt\hbox
to0.0pt{\hss$\scriptstyle\smile$\hss}\kern 2.80048pt}{\scriptscriptstyle\kern
1.42882pt\hbox to0.0pt{\hss$\scriptscriptstyle\mid$\hss}\lower 0.96873pt\hbox
to0.0pt{\hss$\scriptscriptstyle\smile$\hss}\kern 1.42882pt}}^{K}$ be extended
to an automorphism-invariant ternary relation between arbitrary small sets
which satisfies the properties listed in thm:kim-pillay-style?

## References

  * [BY03a] Itay Ben-Yaacov.  Positive model theory and compact abstract theories.  Journal of Mathematical Logic, 03(01):85–118, May 2003.
  * [BY03b] Itay Ben-Yaacov.  Simplicity in compact abstract theories.  Journal of Mathematical Logic, 03(02):163–191, November 2003.
  * [BY03c] Itay Ben-Yaacov.  Thickness, and a categoric view of type-space functors.  Fundamenta Mathematicae, 179:199–224, 2003.
  * [BYP07] Itaï Ben Yaacov and Bruno Poizat.  Fondements de la logique positive.  J. Symbolic Logic 72, pages 1141–1162, 2007.
  * [Cha02] Zoé Chatzidakis.  Properties of forking in $\omega$-free pseudo-algebraically closed fields.  J. Symbolic Logic 67, pages 957––996, 2002.
  * [Cha08] Zoé Chatzidakis.  Independence in (unbounded) PAC fields, and imaginaries.  Around Classification Theory, Leeds, 2008.
  * [CKR20] Artem Chernikov, Byunghan Kim, and Nicholas Ramsey.  Transitivity, lowness, and ranks in NSOP1 theories.  arXiv:2006.10486 [math], June 2020.  arXiv: 2006.10486.
  * [CR16] Artem Chernikov and Nicholas Ramsey.  On model-theoretic tree properties.  Journal of Mathematical Logic 16, May 2016.
  * [d’E21a] Christian d’Elbée.  Forking, imaginaries and other features of ACFG.  Journal of Symbolic Logic, page 1–34, 2021.
  * [d’E21b] Christian d’Elbée.  Generic expansions by a reduct.  Journal of Mathematical Logic, 2021.  published online.
  * [dKN21] Christian d’Elbée, Itay Kaplan, and Leor Neuhauser.  Existentially closed models of fields with a distinguished submodule.  arXiv:2110.02361, 2021.
  * [DKR19] Jan Dobrowolski, Byunghan Kim, and Nicholas Ramsey.  Independence over arbitrary sets in NSOP1 theories.  arXiv:1909.08368 [math], September 2019.
  * [Gra99] Nicholas Granger.  Stability, simplicity, and the model theory of bilinear forms.  PhD thesis, University of Manchester, 1999.
  * [Hay19] Levon Haykazyan.  Spaces of types in positive model theory.  The Journal of Symbolic Logic, 84(2):833–848, June 2019.
  * [HK21] Levon Haykazyan and Jonathan Kirby.  Existentially closed exponential fields.  Israel Journal of Mathematics, 241(1):89–117, March 2021.
  * [Hru] Ehud Hrushovski.  Simplicity and the Lascar group.  Preprint, 1997.
  * [Kam20] Mark Kamsma.  The Kim-Pillay theorem for Abstract Elementary Categories.  The Journal of Symbolic Logic, 85(4):1717–1741, December 2020.
  * [Kim] Byunghan Kim.  Weak canonical bases in NSOP1 theories.  Preprint, 2020.
  * [Kim98] Byunghan Kim.  Forking in Simple Unstable Theories.  Journal of the London Mathematical Society 57, pages 257–267, 1998\.
  * [Kim09] Byunghan Kim.  NTP1 theories.  Slides, Stability Theoretic Methods in Unstable Theories, BIRS, 2009.
  * [KKS14] Byunghan Kim, Hyeungjoon Kim, and Lynn Scow.  Tree indiscernibilities, revisited.  Arch. Math. Logic 53, pages 211–232, 2014.
  * [KP97] Byunghan Kim and Anand Pillay.  Simple theories.  Annals of Pure and Applied Logic 88, pages 149–164, 1997.
  * [KR20] Itay Kaplan and Nicholas Ramsey.  On Kim-independence.  Journal of the European Mathematical Society, 22(5), January 2020\.
  * [KR21] Itay Kaplan and Nicholas Ramsey.  Transitivity of kim-independence.  Advances in Mathematics, 379:107573, 2021.
  * [KRS19] Itay Kaplan, Nicholas Ramsey, and Saharon Shelah.  Local character of kim-independence.  Proceedings of the American Mathematical Society, 147(4):1719–1732, January 2019.
  * [Men20] Rosario Mennuni.  Product of invariant types modulo domination–equivalence.  Arch. Math. Logic 53, pages 1–29, 2020.
  * [Pil00] Anand Pillay.  Forking in the category of existentially closed structures.  Quaderni di Matematica, 6:23–42, 2000.
  * [PY18] Bruno Poizat and Aibat Yeshkeyev.  Positive Jonsson Theories.  Logica Universalis, 12(1):101–127, May 2018.
  * [Sim15] Pierre Simon.  A Guide to NIP Theories.  Lecture Notes in Logic 44. Cambridge University Press, Cambridge, 2015\.


# Polariton-based quantum memristors

Ariel Norambuena  Centro de Investigación DAiTA Lab, Facultad de Estudios
Interdisciplinarios, Universidad Mayor, Chile  e-mail:
<EMAIL_ADDRESS>Felipe Torres  Centro de Nanociencia y
Nanotecnología CEDENNA, Avda. Ecuador 3493, Estación Central, 9170124,
Santiago, Chile  Departamento de Física, Facultad de Ciencias, Universidad de
Chile, Casilla 653, Santiago, 7800024, Chile  Department of Physics,
University of California, San Diego, La Jolla, CA 92093, USA  Massimiliano Di
Ventra  Department of Physics, University of California, San Diego, La Jolla,
CA 92093, USA  Raúl Coto  Centro de Investigación DAiTA Lab, Facultad de
Estudios Interdisciplinarios, Universidad Mayor, Chile  e-mail:
<EMAIL_ADDRESS>
###### Abstract

Information processing and storing by the same physical system is emerging as
a promising alternative to traditional computing platforms. In turn, this
requires the realization of elementary units whose memory content can be
easily tuned and controlled. Here, we introduce a polariton-based quantum
memristor where the memristive nature arises from the inter-cavity polariton
exchange and is controlled by a time-varying atom-cavity detuning. A dynamical
hysteresis is characterized by the fluctuations in the instantaneous polariton
number, where the history information is encoded into a dynamical phase. Using
a Lindblad master equation approach, we find that features of the quantum
memristor dynamics, such as the area and circulation of the hysteresis loop,
showcase a kind of “plasticity” controlled by quantum state initialization.
This makes this quantum memristor very versatile for a wide range of
applications.

## Introduction

Computing in and with memory [1] is emerging as an alternative to our
traditional model of computation, which separates the tasks of storage and
processing of information into two physically distinct units. Memory (time
non-locality) can be realized in many ways and with different types of systems
and devices (see, e.g., Refs. [2] for a review). However, so far particular
attention has been mainly given to classical resistive memories (sometimes
called “memristive elements”) as components of processors with memory [3]. On
the other hand, quantum dynamics may offer additional benefits for the
realization of new storage and information processing capabilities with
applications in quantum information and (neuromorphic) quantum computing [4,
5, 6, 7].

In this respect, initial suggestions of quantum memristive elements have
ranged from superconducting quantum circuits [8, 4, 5, 6] to quantum photonic
devices [9]. In all these cases, the memory mechanism arises from the
combination of quantum feedback [4, 6] and dissipative effects [9]. However,
one can also rely on quasi-particle dynamics, such as low-energy quasi-
particle tunneling [6] and exciton-polariton interactions [10], as other types
of quantum features that can be harnessed to realize quantum memristors.

Among quasi-particles, those originating from strongly correlated light-matter
systems, known as polaritons, offer attractive technological features like
room temperature operation, high dynamical speed, and suitable implementation
[11, 12]. Polaritons have been implemented in a wide variety of systems,
including atoms [13], excitons [14], trapped ions [15] and superconducting
circuits [16], and have allowed the observation of new quantum phases [17, 18,
19, 20].

Recently, more versatile polariton systems have been designed by considering
driving schemes [21, 22, 23]. The external driving sets a different time
scale. Thereby, the system’s time evolution is governed by the difference
between the intrinsic response time scale of the system and the driving time
scale, leading to bistability. Bistability is a hallmark of driven non-linear
classical and quantum systems [24, 25, 26, 27, 28, 29]. In the latter, this is
due to the presence of two stable branches where quantum fluctuations switch
the states of the system between these two stable solutions [29, 30].
Furthermore, the endeavors to address the response of a system under the
influence of an external driving focusing on understanding the out-of-
equilibrium dynamics have opened new perspectives. In such direction, we found
the control of Mott-insulator and superfluid states [31] and the appearance of
dynamical hysteresis in photonic devices [29, 32].

In this work, we introduce a polariton-based quantum memristor (PQM) whose
behavior is controlled by a time-dependent atomic modulation providing a
hysteretic response over the variance of the number of quasi-particles. The
underlying mechanism relies on the inter-cavity polariton exchange (hopping)
and atomic modulation. We demonstrate that memory effects in our PQM depend on
the initialization of the system (plasticity) and the time scale of the
driving. This minimalistic setup exploits the inter-cavity polariton exchange,
leading to interesting phenomena in more complex optical arrays.

Figure 1: Schematic representation of the polariton-based quantum memristor
composed of two coupled cavities QED with atomic modulation $\xi(t)$. Each
cavity interacts with a two-level atom, and photon hopping is allowed between
the cavities with coupling strength $J$. Furthermore, each cavity experiences
atomic $\gamma_{i}^{\rm a}$ and photonic $\gamma_{i}^{\rm c}$ losses.

## Polariton-based quantum memristor

The envisioned structure consists of two interacting optical cavities driven
by a time-dependent atomic modulation $\xi(t)$. Figure 1 shows a schematic
representation of the PQM composed of two cavities doped with a two-level
system (we label $\varphi_{g}$ and $\varphi_{e}$), where the intra-cavity
light-matter coupling allows the generation of cavity QED polaritons quasi-
particles, and the photon hopping between cavities leads to polariton
exchange. Dissipative effects arise from imperfect mirrors with photonic decay
$(\gamma_{i}^{\rm c})$ and spontaneous emission with atomic decay
$(\gamma_{i}^{\rm a})$ per site ($i=1,2$). The modulation serves different
purposes such as manipulation of qubits [31, 33], sideband transitions [34]
and control of atomic quantum gases [35]. The dynamics of the light-matter
interaction is given by the Jaynes-Cummings-Hubbard (JCH) Hamiltonian,

$H_{0}=\sum_{i=1}^{2}\left[\omega^{\rm c}a_{i}^{\dagger}a_{i}+\omega^{\rm
a}\sigma_{i}^{\dagger}\sigma_{i}+g\left(a_{i}\sigma_{i}^{\dagger}+a_{i}^{\dagger}\sigma_{i}\right)\right]-J\left(a_{1}a_{2}^{\dagger}+a_{2}^{\dagger}a_{1}\right),$
(1)

here $\sigma_{i}=|\varphi_{e}\rangle_{i}\langle\varphi_{g}|_{i}$ are Pauli
operators for the two-level atom at the $i$-th cavity, $a_{i}^{\dagger}$ and
$a_{i}$ are photon creation and annihilation operators, respectively,
$\omega^{\rm c}$ ($\omega^{\rm a}$) is the cavity (atom) frequency, $g$ is the
atom-field coupling strength, and $J$ is the photon hopping strength between
neighboring cavities. In this work, we set $\omega^{\rm c}=\omega^{\rm a}$
unless stated otherwise. The total Hamiltonian of the system is modelled as

$H(t)=H_{0}+\xi(t)\sum^{2}_{i=1}\sigma_{i}^{\dagger}\sigma_{i}.$ (2)

We investigate the dynamics induced by a Gaussian modulation given by

$\xi(t)=\xi_{i}+(\xi_{f}-\xi_{i})e^{-(t-T)^{2}/2\sigma_{w}^{2}},$ (3)

where $\xi_{f}>\xi_{i}>0$ and $\sigma_{w}$ is the characteristic width of the
Gaussian modulation with the property. We explore the dynamics in the regime
$T^{2}/(2\sigma_{w}^{2})>1$ leading to a Gaussian profile satisfying
$\xi(0)=\xi(2T)\approx\xi_{i}$. We introduce the instantaneous polaritonic
basis
$|n+\rangle=\sin\theta_{n}(t)|n,\varphi_{g}\rangle+\cos\theta_{n}(t)|n-1,\varphi_{e}\rangle$
and
$|n-\rangle=\cos\theta_{n}(t)|n,\varphi_{g}\rangle-\sin\theta_{n}(t)|n-1,\varphi_{e}\rangle$
with $\theta_{n}(t)=(1/2)\mbox{tan}^{-1}(2\sqrt{n}g/\xi(t))$, and
instantaneous eigenenergies $E_{n\pm}(t)=\omega^{\rm
c}n+\xi(t)/2\pm(\xi^{2}(t)+4g^{2}n)^{1/2}/2$, where $n$ is the number of
photons per cavity. Moreover, sudden inter-branch transitions could be
activated if $\xi(t)$ has some non-adiabatic features (see Supplementary
Information for further details).

Figure 2: (a) Dynamical hysteresis curves for a closed PQM with $g=1$,
$J=10^{-2}g$, $\xi_{i}=10g$, $\xi_{f}=10^{3}g$, $\sigma_{w}=T/4$, $\omega^{\rm
c}=10^{4}g$, and three different times $T$ in the Gaussian modulation. The
initial state is a Mott-insulator-like state
$|\Psi(0)\rangle=|1-\rangle\otimes|1-\rangle$ for $\xi(t=0)=\xi_{i}$. (b)
Population dynamics. Solid and dashed lines correspond to the exact numerical
solutions for the Gaussian time-dependent and a large constant detunings,
respectively. The case $\xi(t)=\xi_{f}\gg g$ shows differences with the time-
dependent case, since dynamical phases introduce a small shift in the
populations.

When the atomic and photonic frequencies are similar this system exhibit a
Mott-phase characterized by localized excitations caused by the photon
blockade, i.e., the effective repulsion between the excitations [36, 13], that
prevents fluctuations of the total number of polaritonic excitations
$\sum_{i}N_{i}$ with
$N_{i}=a_{i}^{\dagger}a_{i}+\sigma_{i}^{\dagger}\sigma_{i}$ the number of
quantum excitations in the $i$th cavity. A photonic superfluid phase,
characterized by delocalized excitations, emerges from the detuning between
atomic and photonic frequencies resembling the flow of charge in a classical
resistive memory. During the Mott-superfluid phase transition, the
delocalization of the excitations and the polariton exchange lead to an
effective quasi-particle flow throughout the optical cavities. In our system,
time non-locality arises from the atomic modulation, and memory effects are
studied using the variance in the polariton number
$\text{var}[N_{i}(t)]=\mbox{Tr}\left[N_{i}^{2}\rho(t)\right]-\mbox{Tr}\left[N_{i}\rho(t)\right]^{2}$.

We describe the open dynamics using the Markovian master equation
$\dot{\rho}=-i[H(t),\rho]+\mathcal{D}\rho$, where losses are modeled using the
Lindblad form:

$\displaystyle{\cal D}\rho$ $\displaystyle=$
$\displaystyle\sum^{2}_{i=1}\left[\gamma_{i}^{\rm c}a_{i}\rho
a_{i}^{\dagger}+\gamma_{i}^{\rm
a}\sigma_{i}\rho\sigma_{i}^{\dagger}-\frac{1}{2}\left\\{a_{i}^{\dagger}a_{i}+\sigma_{i}^{\dagger}\sigma_{i},\rho\right\\}\right].$

For concreteness, we first consider an initial Mott-insulator-like state
$|\Psi(0)\rangle=|1-\rangle^{\xi(0)}\otimes|1-\rangle^{\xi(0)}$, where
$|1-\rangle^{\xi(0)}$ means that the initial state is calculated for $\xi(0)$.
Additionally, we set $\xi(t)\geq 10g$ in Eq. (3), which serves two purposes.
Firstly, it enables inter-cavity polariton exchange by suppressing the photon-
blockade. Secondly, variations on $\theta_{n}(t)$ are small, constraining the
effect of the modulation to small shifts on the variance. In Fig. 2(a) we show
that the evolution of $\text{var}[N_{i}(t)]$ describes a hysteresis loop as a
function of $\xi(t)$. We remark that the transition Mott-Superfluid takes
place regardless the instantaneous value of $\xi(t)$. However, the driving
time scale originates a shift on the population, which can be observed in Fig
2(b) and Eqs.(S11)-(S12) in the Supplementary Information. Therefore, as the
position of the Gaussian peak ($T$) changes different loops appear, but there
is a critical value ($\tau=\pi/(4J)\approx 78.5g^{-1}$ for $J=10^{-2}g$) that
suppress the hysteretic behavior (see Supplementary Information for the
estimation of the critical time). To understand this, we provide a simple
argument based on the time-symmetry of the functions involved.

Let us consider the time evolution of the occupation probabilities depicted in
Fig. 2(b). Throughout the manuscript, we refer to the occupation probability
for each instantaneous two-body polaritonic state $|n\pm,m\pm\rangle$ as
$p_{n\pm,m\pm}=\langle n\pm,m\pm|\rho(t)|n\pm,m\pm\rangle$. It is worth
noticing that there is a characteristic time ($\tau$) that corresponds to the
condition $p_{1-,1-}(\tau)=0$, where the system is in a superfluid-like state,
i.e. only states $|2-,0g\rangle$ and $|2-,0g\rangle$ (with the same
probabilities) are present. Around this time ($\tau$) the dynamics is roughly
time-symmetric, $\text{var}[N_{i}(\tau+t)]\approx\text{var}[N_{i}(\tau-t)]$.
Moreover, the detuning in Eq. (3) is a Gaussian function, and it is symmetric
around $T$. Hence, when $T=\tau$, $\text{var}[N_{i}(t)]$ and $\xi(t)$ exhibit
the same time-symmetry which makes the hysteresis loop vanish, see Fig. 2(a).
However, as $\tau$ and $T$ differ, the hysteresis loop is restored. The
hysteresis loop arises from the delay between the driving time ($T$) and the
time-scale of the system ($\tau$). This implies that the irreversibility of
the dynamical process is produced by the response time scale of the system
against external driving. We analyze the hysteresis loop area and its
circulation for different initial conditions in the closed and open dynamics.
In the Supplementary Information we analyze the characteristic time $\tau$ in
an effective spin-$1$ system that mimics our PQM.

Figure 3: Hysteresis curves for the Mott insulator (superfluid) initial
condition $|\Psi(0)\rangle=|1-,1-\rangle$
($|\Psi(0)\rangle=(1/\sqrt{2})[|2-,0g\rangle+|0g,2-\rangle]$) for the closed
((a) and (c)) and open ((b) and (d)) dynamics of the PQM. For all simulations
the parameters are $g=1$, $J=10^{-2}g$, $\xi_{i}=10g$, $\xi_{f}=10^{3}g$,
$\sigma_{w}=T/4$, and $\omega^{\rm c}=10^{4}g$. For the open case we set
$\gamma_{i}^{\rm a,c}=10^{-2}g$ for each cavity $i=1,2$.

## Hysteresis loop area and circulation

The oriented area of a simple connected curve $(x(t),y(t))$ can be defined as
$-\oint ydx$. In our particular case, we have $x(t)=\xi(t)$ (control) and
$y(t)=\mbox{var}[N_{i}]$ (response). The enclosed area in the absence of
dissipation is approximately given by

$A\approx{\eta\over\sqrt{2}}(\xi_{f}-\xi_{i})\pi^{3/2}{\sigma_{w}\over\tau}\sin\left({\pi
T\over\tau}\right)e^{-(\pi\sigma_{w}/(\sqrt{2}\tau))^{2}},$ (5)

where $\eta=1-2\mbox{var}[N_{i}(0)]$ is a factor that depends on the state
initialization (see Supplementary Information). We first notice that $A=0$ for
$T=\tau$ ($\sin(\pi T/\tau)=0$), which explains why the hysteresis vanishes in
Fig. 2(a) dotted line. When $T\neq\tau$, the sign of $A$ accounts for the
direction of the circulation, where $A>0$ (anti-clockwise) and $A<0$
(clockwise). We analytically find that the sign of the area depends on three
factors, namely, the variation $\xi_{f}-\xi_{i}$ (which is positive in our
work), the ratio $T/\tau$, and the factor $\eta$ that accounts for state
initialization. We show in Fig. 3 that the inter-cavity polariton exchange and
the time-dependent driving $\xi(t)$ induce different responses in the variance
of the polariton number $\text{var}[N_{i}]$ depending on the initial state
$\rho(0)=|\Psi(0)\rangle\langle\Psi(0)|$.

In Fig. 3(a) we show the case where $\eta=1$ ($\mbox{var}[N_{i}(0)]=0$),
leading to anti-clockwise (clockwise) circulation for $T<\tau$ ($T>\tau$). In
Fig. 3(c), instead, we show the case $\eta=-1$ ($\mbox{var}[N_{i}(0)]=1$), for
which the circulation has been inverted with respect to the previous case. In
Fig. 3(b)-(d) we show the hysteresis for the open dynamics. The main
difference is that there is still a loop in Fig. 3(b), while in the Fig. 3(d),
the loop breaks. The latter is a consequence of the initial condition that
sets the initial variance ($\text{var}[N_{i}(0)]$) at the maximum value, while
losses do not allow the system to reach that state again. Finally, we remark
that the dynamics is no longer symmetric around $\tau=78.5g^{-1}$, thus the
condition $T=\tau$ for shrinking the area fails.

Figure 4: (a) Hysteresis curves for the Jaynes-Cummings modulated with a time-
dependent laser with frequency $\omega_{l}$ and a Gaussian laser amplitude
$I(t)$. The hysteresis curves are plotted for the mean number of photons
$\langle a^{\dagger}a\rangle$ as a function of the laser amplitude $I(t)/g$
using different polaritonic initial conditions $|\Psi(0)\rangle$. For the
simulation we use the parameters $g=1$, $w^{\rm a}=w^{\rm c}=w^{\rm
l}=10^{4}g$, $\gamma^{\rm a}=\gamma^{\rm c}=g/10$ (photonic and atomic losses)
and $T=0.95\tau$. (b) Hysteresis curves for the non-linear Kerr model
introduced in Ref. [32] with time-dependent laser amplitude $F(t)$ following a
triangular shape with sweep time $t_{s}=10\gamma^{-2}\Delta F$, $\Delta
F=3F_{0}$, $F_{0}=\gamma$. The mean number of photons $\langle
a^{\dagger}a\rangle$ as a function of $F/\gamma$ shows a closed hysteresis
loop. For the simulations we consider $\gamma=1$, $\Delta=2\gamma$,
$U=0.1\gamma$, $\gamma=1$ and different initial Fock states $|\Psi(0)\rangle$.
(c) Periodic modulation $\xi(t)$ and its effects on the hysteresis curve
$(\xi(t)/g,\mbox{var}[N_{i}])$ as a function of time $t/T$. The 3D plot show
how the response function $\mbox{var}[N_{i}]$ describes a spiral when the
train-pulse $\xi(t)$ is introduced. For the simulation we consider a closed
dynamics with $g=1$, $J=10^{-2}g$, $\omega^{\rm c}=10^{4}g$, $T=0.95\tau$,
$\sigma_{w}=T/4$, and the initial condition $|\Psi(0)\rangle=|1-,1-\rangle$.
(d) PQM dynamics simulated from Eq. (9) starting from the initial conditions
$|\Psi(0)\rangle=|1-,1-\rangle$ (solid line) and
$|\Psi(0)\rangle=(|2-,0g\rangle+|0g,2-\rangle)/\sqrt{2}$. The inset shows the
PQM circuit analogy, Eqs (9) and (10). For simulations we use the parameters
$g=1$, $J=10^{-2}g$, $\xi_{i}=10g$, $\xi_{f}=10^{3}g$, $\sigma=T/4$, and
$\omega^{\rm c}=10^{4}g$, and $\gamma_{i}^{\rm a,c}=0$ for each cavity
$i=1,2$.

Quantum platforms as non-equilibrium photonic systems also display hysteresis
loops as a function of the photonic number [32]. However, these systems are
independent of the initial condition. For instance, consider a single cavity
driven by a continuous wave with frequency $\omega^{p}$ and coupling strength
$I(t)$, while the atom is driven with frequency $\omega^{l}$ and Rabi coupling
$\Omega$. In a multi-rotating frame with the atom and cavity frequencies, we
get

$\displaystyle H_{\rm
1}(t)=\Delta_{a}\sigma^{\dagger}\sigma+\Delta_{c}a^{\dagger}a+g\left(a^{\dagger}\sigma
e^{i\Delta_{1}t}+a\sigma^{\dagger}e^{-i\Delta_{1}t}\right)+i\Omega(\sigma^{\dagger}-\sigma)+iI(t)\left(a^{\dagger}-a\right),$
(6)

where $\Delta_{a}=\omega^{\rm a}-\omega^{\rm l}$, $\Delta_{c}=\omega^{\rm
c}-\omega^{\rm p}$, and $\Delta_{1}=\omega^{\rm p}-\omega^{\rm l}$. Photonic
and atomic losses for individual cavities are described by the operator
$\mathcal{D}\rho$ defined in Eq. (Polariton-based quantum memristor). For
simplicity we set $\Omega=10^{-6}g\ll 1$, and we model
$I(t)=I_{0}e^{-(t-T)^{2}/2\sigma_{w}^{2}}$ with a Gaussian profile by setting
$I_{0}=g$, $\sigma=T/4$, and $T=10^{3}g^{-1}$. In Fig. 4 (a) we show that the
system undergoes a hysteretic behavior for the resonant condition
$\Delta_{a}=\Delta_{c}=\Delta_{1}=0$, but the loops remain unaffected by the
initial conditions. In addition, we consider the driven non-linear Kerr model,
described in Ref [32], where the system Hamiltonian in the rotating frame at
the pumping frequency $\omega^{\rm p}$ is:

$H_{\rm KM}(t)=-\Delta
a^{\dagger}a+\frac{U}{2}(a^{\dagger})^{2}a^{2}+F(t)(a^{\dagger}+a),$ (7)

where $\Delta=\omega^{\rm p}-\omega^{\rm a}$ and $F(t)=F_{0}+(t/t_{s})\Delta
F\theta(t_{s}-t)-(t-2t_{s})/t_{s}\Delta F\theta(t-t_{s})$ represent a
triangular shape, with $\Delta F$ the maximum and $\theta(t)$ the Heaviside
step function. For simulating the non-linear Kerr model we use the master
equation $\dot{\rho}=-i[H_{\rm KM}(t),\rho]+\gamma/2(2a\rho
a^{\dagger}-\\{a^{\dagger}a,\rho\\})$. Once again, we observe dynamical
hysteresis without a dependency on the initial conditions, see Fig. 4 (b). To
illustrate our results further, we consider a pulse train for the external
modulation $\xi(t)$. After each pulse, the system accumulates a phase due to
the shifts in the population (see Fig. 2(b)), and ends in a different state,
which serves as the initial condition for the incoming pulse. Henceforth, we
find in Fig. 4 (c) a spiral time-dependent hysteresis loops. We use the term
plasticity to denote this particular dynamical sensitiveness to initial states
in our scheme.

Having demonstrated the notion of plasticity, a natural question arises: can
our model be mathematically connected with the underlying physics of circuit
elements with memory? In particular, this will be useful for physical
implementations of quantum neuromorphic computing [4, 7] which simulate
artificial or biological neural networks.

## Circuit analogy

In classical resistive memories [2], the current-voltage (I-V) relation
$I(t)=G(x(t),V(t),t)V(t)$ holds, where $G$ is the memory conductance, and the
internal state variables $x(t)$ track the past state history of the system
according to the applied input signal. The inclusion of an updating function
$\dot{x}(t)=f(x(t),V(t),t)$ typically yields non-linear dynamics under
periodic driving and produces an I-V hysteresis loop [2]. To make a connection
between the above circuit elements and our PQM, we use the time-dependent
variation of the polariton variance and the Lindblad master equation to obtain

$\frac{dy}{d\xi}\dot{\xi}+b(\xi,\rho,t)\xi(t)=a(t,\xi,\rho,t),$ (8)

where $y(t)=\text{Re}\\{\text{var}[N_{i}]\\}$ accounts for the real part of
the response function,
$a(t,\xi,\rho,t)=\text{Re}\\{-i\langle[\alpha_{i},H_{0}]\rangle+\text{Tr}[\alpha_{i}{\cal
D}\rho]\\}$, $b(t)=\text{Re}\\{i\langle[\alpha_{i},H_{i}]\rangle\\}$ are
dynamical functions that depend on the input function $\xi(t)$, the state of
the system $\rho(t)$ and the Hamiltonians $H_{0}$ and
$H_{i}=\sum_{i}\sigma_{i}^{\dagger}\sigma_{i}$. Here,
$\alpha_{i}=(N_{i}-\langle N_{i}\rangle)^{2}-\langle N_{i}\rangle^{2}$ with
$\langle N_{i}\rangle=\text{Tr}(N_{i}\rho)$. We note that the form of the
dynamical system given in Eq. (8) resembles a forced RC circuit, where the
following analogies can be recognized: $i(t)\rightarrow\dot{\xi}$ (current),
$R\rightarrow dy/d\xi$ (resistance), $C\rightarrow b$ (capacitance), and
$\varepsilon\rightarrow a$ (voltage), see Fig. 4(d). If we define
$R(\xi,\rho,t)=dy/d\xi$ as differential resistance, we obtain the following
Ohmic-like law:

$\displaystyle V_{R}$ $\displaystyle=$ $\displaystyle R(\xi,\rho,t)\dot{\xi},$
(9) $\displaystyle\dot{\rho}$ $\displaystyle=$ $\displaystyle f(\rho,\xi,t).$
(10)

As shown in Fig. 4 (d), the above system mimics the dynamics of a current-
controlled memristive system [4], where the voltage-like $V_{R}=\dot{y}$
varies as a function of the current $\dot{\xi}$, and the internal state
evolves following a Markovian model $f(\rho,\xi,t)=-i[H(t),\rho]+{\cal
D}\rho$. We remark that $R(\xi,\rho,t)=dy/d\xi$ can be positive, zero or even
negative, and such behavior is connected with the initial condition and the
relative time $\tau-T$. For instance, as shown in Fig. 3 (a) for $T=0.95\tau$,
the variance starts from $y(t=0)=0$ and increases as a function of the
detuning $\xi$, which implies that initially $R(\xi,\rho,t)=dy/d\xi>0$. On the
contrary, in Fig. 3 (c) for $T=0.95\tau$, the initial value is $y(t=0)=1$ and
the variance decreases in terms of $\xi$, thus $R(\xi,\rho,t)=dy/d\xi<0$.

Table 1: Characteristics of different types of quantum memristors. Quantum
Memristors
---
System |  Mechanisms |  Features |  Implementation Status
Cavity QED polariton |  polariton exchange and external driving |  Sensitiveness to initial states, controlled memory effect without training, robust behavior against dissipative effects . |  Theoretical
Beam splitter with a tunable reflectivity[37] |  Classical feedback |  Resistive Memory, fast update speed |  Experimental
Superconducting Circuits |  Quantum feedback [4] or quasi-particle tunneling [6] |  Sensitiveness to system initialization [6], Low temperature management, training, |  Theoretical
Double quantum dot[38] |  Capacitive coupling |  Low temperature management, small tunnel barrier |  Theoretical

As an example, let us consider two initial conditions
$|\Psi(0)\rangle=|1-,1-\rangle$ (solid line with resistance $R_{1}$) and
$|\Psi(0)\rangle=(|2-,0g\rangle+|0g,2-\rangle)/\sqrt{2}$ (dashed line with
resistance $R_{2}$). In both cases, we have $V_{R}(t=0)=0$ and
$\gamma_{i}^{\rm a,c}=0$ ($i=1,2$). The simulations illustrated in Fig. 4 (d)
show two pinched hysteresis loops for the voltage $V_{R}$ as a function of the
current $\dot{\xi}$, where both curves phenomenologically describe the
memristive element physics. Finally, we numerically corroborate that the only
difference between the two curves in Fig. 4 (d) is the sign of the
resistances, i.e. $R_{1}=-R_{2}$, which explains the asymmetrical aspect.
Finally, there have been a few suggestions of quantum memristors. In Table1 we
report a comparative analysis between currently proposed implementations. The
PQM presented here shares some features with superconducting circuits, in the
sense that both use quasi-particle manipulation while offering remarkable
tunability.

## Conclusions

We have proposed a novel polariton-based quantum memristor (PQM), where the
memristive behavior arises from the inter-cavity polariton exchange and is
controlled by the time-varying atom-cavity detuning (atomic modulation). The
hysteresis loop area and circulation depend on the quantum state
initialization and time scales between the system and modulation. This
behavior emulates a sort of “plasticity of memristive devices”. Remarkably,
this plasticity is not observed in hysteretic systems such as the single
Jaynes-Cummings and non-linear Kerr model. This suggests that quasi-particle
dynamics offer additional functionalities. From a technical point of view,
implementing an elementary cavity QED polariton quantum memristor, whose
memory content is controlled by the state initialization, adds interesting
features like room temperature operation and high dynamical speed. However,
perhaps the most significant aspect of PQMs, as opposed to classical
memristive elements, may be the energy-efficient flow of quasi-particles that
does not involve thermal dissipation due to the Joule heating. This essential
feature is crucial for low-energy applications, and it makes the PQM an
exciting and versatile candidate for the use of strongly-correlated light-
matter systems in the field of memcomputing and (neuromorphic) quantum
computation.

## References

  * [1] Di Ventra, M. & Pershin, Y. V.  The parallel approach.  _Nature Physics_ 9, 200 (2013).
  * [2] Pershin, Y. V. & Di Ventra, M.  Memory effects in complex materials and nanoscale systems.  _Advances in Physics_ 60, 145–227, DOI: 10.1080/00018732.2010.544961 (2011).
  * [3] Traversa, F. L. & Di Ventra, M.  Universal MemComputing machines.  _IEEE Trans. Neural Netw. Learn. Syst._ 26, 2702, DOI: 10.1109/TNNLS.2015.2391182 (2015).
  * [4] Pfeiffer, P., Egusquiza, I. L., Di Ventra, M., Sanz, M. & Solano, E.  Quantum memristors.  _Scientific Reports_ 6, 29507, DOI: 10.1038/srep29507 (2016).
  * [5] Gonzalez-Raya, T., C., Lukens, J. M. & Sanz, M.  Quantum memristors in frequency-entangled optical fields.  _Materials_ 13 (2020).
  * [6] Salmilehto, J., Deppe, F., Di Ventra, M., Sanz, M. and Solano, E.  Quantum memristors with superconducting circuits.  _Scientific Reports_ 7, 42044, DOI: 10.1038/srep42044 (2017).
  * [7] Marković, D. & Grollier, J.  Quantum neuromorphic computing.  _Applied Physics Letters_ 117, 150501, DOI: 10.1063/5.0020014 (2020).  https://doi.org/10.1063/5.0020014.
  * [8] Peotta, S. & Di Ventra, M.  Superconducting memristors.  _Physical Review Applied_ 2, 034011 (2014).
  * [9] M. Sanz, E. S., L. Lamata.  Invited article: Quantum memristors in quantum photonics.  _APL Photonics_ 3, 080801, DOI: 10.1063/1.5036596 (2018).
  * [10] Mirek, R., P., Matuszewski, M. & et. al.  Neuromorphic binarized polariton networks.  _Nano Letters_ 21, 3715–3720, DOI: 10.1021/acs.nanolett.0c04696 (2021).
  * [11] Boulier, T. _et al._ Microcavity polaritons for quantum simulation.  _Advanced Quantum Technologies_ 2000052, DOI: 10.1002/qute.202000052 (2020).  https://onlinelibrary.wiley.com/doi/pdf/10.1002/qute.202000052.
  * [12] Basov, D. N., Asenjo-Garcia, A., Schuck, P. J., Zhu, X. & Rubio, A.  Polariton panorama:.  _Nanophotonics_ 10, 549–577, DOI: doi:10.1515/nanoph-2020-0449 (2021).
  * [13] Birnbaum, K. M. _et al._ Photon blockade in an optical cavity with one trapped atom.  _Nature_ 436, 87–90, DOI: 10.1038/nature03804 (2005).
  * [14] Weisbuch, C., Nishioka, M., Ishikawa, A. & Arakawa, Y.  Observation of the coupled exciton-photon mode splitting in a semiconductor quantum microcavity.  _Phys. Rev. Lett._ 69, 3314–3317, DOI: 10.1103/PhysRevLett.69.3314 (1992).
  * [15] Toyoda, K., Matsuno, Y., Noguchi, A., Haze, S. & Urabe, S.  Experimental realization of a quantum phase transition of polaritonic excitations.  _Phys. Rev. Lett._ 111, 160501, DOI: 10.1103/PhysRevLett.111.160501 (2013).
  * [16] Raftery, J., Sadri, D., Schmidt, S., Türeci, H. E. & Houck, A. A.  Observation of a dissipation-induced classical to quantum transition.  _Phys. Rev. X_ 4, 031043, DOI: 10.1103/PhysRevX.4.031043 (2014).
  * [17] Hartmann, M. J., Brandão, F. G. S. L. & Plenio, M. B.  Strongly interacting polaritons in coupled arrays of cavities.  _Nature Physics_ 2, 849–855, DOI: 10.1038/nphys462 (2006).
  * [18] Greentree, A. D., Tahan, C., Cole, J. H. & Hollenberg, L. C. L.  Quantum phase transitions of light.  _Nature Physics_ 2, DOI: 10.1038/nphys466 (2006).
  * [19] Fitzpatrick, M., Sundaresan, N. M., Li, A. C. Y., Koch, J. & Houck, A. A.  Observation of a dissipative phase transition in a one-dimensional circuit qed lattice.  _Phys. Rev. X_ 7, 011016, DOI: 10.1103/PhysRevX.7.011016 (2017).
  * [20] Kasprzak, J. _et al._ Bose–einstein condensation of exciton polaritons.  _Nature_ 443, 409–414, DOI: 10.1038/nature05131 (2006).
  * [21] Johansen, C. H., Lang, J. & Piazza, F.  Multimode-floquet-polariton superradiance (2020).  2011.12309.
  * [22] Clark, L. W. _et al._ Interacting floquet polaritons.  _Nature_ 571, 532–536, DOI: 10.1038/s41586-019-1354-5 (2019).
  * [23] Mundada, P. S., Gyenis, A., Huang, Z., Koch, J. & Houck, A. A.  Floquet-engineered enhancement of coherence times in a driven fluxonium qubit.  _Phys. Rev. Applied_ 14, 054033, DOI: 10.1103/PhysRevApplied.14.054033 (2020).
  * [24] Gibbs, H. M., McCall, S. L. & Venkatesan, T. N. C.  Differential gain and bistability using a sodium-filled fabry-perot interferometer.  _Phys. Rev. Lett._ 36, 1135–1138, DOI: 10.1103/PhysRevLett.36.1135 (1976).
  * [25] Felber, F. S. & Marburger, J. H.  Theory of nonresonant multistable optical devices.  _Applied Physics Letters_ 28, 731–733, DOI: 10.1063/1.88632 (1976).  https://doi.org/10.1063/1.88632.
  * [26] Lugiato, L. A.  Ii theory of optical bistability.  vol. 21 of _Progress in Optics_ , 69–216, DOI: https://doi.org/10.1016/S0079-6638(08)70122-7 (Elsevier, 1984).
  * [27] Rempe, G., Thompson, R. J., Brecha, R. J., Lee, W. D. & Kimble, H. J.  Optical bistability and photon statistics in cavity quantum electrodynamics.  _Phys. Rev. Lett._ 67, 1727–1730, DOI: 10.1103/PhysRevLett.67.1727 (1991).
  * [28] Gripp, J., Mielke, S. L., Orozco, L. A. & Carmichael, H. J.  Anharmonicity of the vacuum rabi peaks in a many-atom system.  _Phys. Rev. A_ 54, R3746–R3749, DOI: 10.1103/PhysRevA.54.R3746 (1996).
  * [29] Rodriguez, S. R. K. _et al._ Probing a dissipative phase transition via dynamical optical hysteresis.  _Phys. Rev. Lett._ 118, 247402, DOI: 10.1103/PhysRevLett.118.247402 (2017).
  * [30] Huybrechts, D. & Wouters, M.  Dynamical hysteresis properties of the driven-dissipative bose-hubbard model with a gutzwiller monte carlo approach.  _Phys. Rev. A_ 102, 053706, DOI: 10.1103/PhysRevA.102.053706 (2020).
  * [31] Tancara, D. _et al._ Steering interchange of polariton branches via coherent and incoherent dynamics.  _Phys. Rev. A_ 103, 053708, DOI: 10.1103/PhysRevA.103.053708 (2021).
  * [32] Casteels, W., Storme, F., Le Boité, A. & Ciuti, C.  Power laws in the dynamic hysteresis of quantum nonlinear photonic resonators.  _Phys. Rev. A_ 93, 033824, DOI: 10.1103/PhysRevA.93.033824 (2016).
  * [33] Strand, J. D. _et al._ First-order sideband transitions with flux-driven asymmetric transmon qubits.  _Phys. Rev. B_ 87, 220505, DOI: 10.1103/PhysRevB.87.220505 (2013).
  * [34] Beaudoin, F., da Silva, M. P., Dutton, Z. & Blais, A.  First-order sidebands in circuit qed using qubit frequency modulation.  _Phys. Rev. A_ 86, 022305, DOI: 10.1103/PhysRevA.86.022305 (2012).
  * [35] Eckardt, A.  Colloquium: Atomic quantum gases in periodically driven optical lattices.  _Rev. Mod. Phys._ 89, 011004, DOI: 10.1103/RevModPhys.89.011004 (2017).
  * [36] Angelakis, D. G., Santos, M. F. & Bose, S.  Photon-blockade-induced mott transitions and $xy$ spin models in coupled cavity arrays.  _Phys. Rev. A_ 76, 031805, DOI: 10.1103/PhysRevA.76.031805 (2007).
  * [37] Gonzalez-Raya, T., Lukens, J. M., Céleri, L. C. & Sanz, M.  Quantum memristors in frequency-entangled optical fields.  _Materials_ 13, DOI: 10.3390/ma13040864 (2020).
  * [38] Li, Y. _et al._ Double quantum dot memristor.  _Phys. Rev. B_ 96, 075446, DOI: 10.1103/PhysRevB.96.075446 (2017).

## Acknowledgements

The authors thank Peter Rabl for helpful insights and comments on the
manuscript. A.N. acknowledges financial support from Universidad Mayor through
the Postdoctoral Fellowship. RC acknowledges financial support from Fondecyt
Iniciación No. 11180143. F. T. acknowledges financial support from Fondo
Nacional de Investigaciones Científicas y Tecnológicas (FONDECYT, Chile) under
grants 1211902 and Centro de Nanociencia y Nanoteccnología CEDENNA,
Financiamiento Basal para Centros Científicos y Tecnológicos de Excelencia
AFB180001. M.D. acknowledges financial support from the Department of Energy
under Grant No. DE-SC0020892.

## Author contributions

A. N. and R. C. created the theoretical model. F. T. conceived the connection
with memristive elements, and M. D. realized insightful comments and
observations related to memristor physics. A. N. performed all numerical
simulations. All authors prepared the manuscript, proofread the paper, made
comments, and approved the manuscript.

## Competing interests

M.D. is the co-founder of MemComputing, Inc. (https://memcpu.com/) that is
attempting to commercialize the memcomputing technology. All other authors
declare no competing interests.