# Human Guided Ground-truth Generation for Realistic Image Super-resolution

Du Chen1, Jie Liang1,2∗, Xindong Zhang1,2, Ming Liu1,3, Hui Zeng2 and Lei
Zhang1,2
1The Hong Kong Polytechnic University, 2OPPO Research Institute, 3Harbin
Institute of Technology
{csdud.chen<EMAIL_ADDRESS>{liang27jie<EMAIL_ADDRESS>
{csxdzhang<EMAIL_ADDRESS>Equal contribution.Corresponding
author. This work is supported by the Hong Kong RGC RIF grant (R5001-18) and
the PolyU-OPPO Joint Innovation Lab.

###### Abstract

How to generate the ground-truth (GT) image is a critical issue for training
realistic image super-resolution (Real-ISR) models. Existing methods mostly
take a set of high-resolution (HR) images as GTs and apply various
degradations to simulate their low-resolution (LR) counterparts. Though great
progress has been achieved, such an LR-HR pair generation scheme has several
limitations. First, the perceptual quality of HR images may not be high
enough, limiting the quality of Real-ISR outputs. Second, existing schemes do
not consider much human perception in GT generation, and the trained models
tend to produce over-smoothed results or unpleasant artifacts. With the above
considerations, we propose a human guided GT generation scheme. We first
elaborately train multiple image enhancement models to improve the perceptual
quality of HR images, and enable one LR image having multiple HR counterparts.
Human subjects are then involved to annotate the high quality regions among
the enhanced HR images as GTs, and label the regions with unpleasant artifacts
as negative samples. A human guided GT image dataset with both positive and
negative samples is then constructed, and a loss function is proposed to train
the Real-ISR models. Experiments show that the Real-ISR models trained on our
dataset can produce perceptually more realistic results with less artifacts.
Dataset and codes can be found at https://github.com/ChrisDud0257/HGGT.

##  1 Introduction

Figure 1: From left to right and top to bottom: one original HR image (Ori) in
the DIV2K[1] dataset, two of its enhanced positive versions (Pos-1 and Pos-2)
and one negative version (Neg). The positive versions generally have clearer
details and better perceptual quality, while the negative version has some
unpleasant visual artifacts. Please zoom in for better observation.

Owing to the rapid development of deep learning techniques[14, 22, 44, 19,
18], the recent years have witnessed the great progress in image super-
resolution (ISR)[8, 9, 10, 23, 26, 33, 27, 56, 52, 28, 48, 12, 2, 13, 35, 51,
46, 29, 45, 31, 32, 54], which aims at generating a high-resolution (HR)
version of the low-resolution (LR) input. Most of the ISR models (_e.g_.,
CNN[38, 37] or transformer[5, 29] based ones) are trained on a large amount of
LR-HR image pairs, while the generation of LR-HR image pairs is critical to
the real-world performance of ISR models.

Most of the existing ISR methods take the HR images (or after some sharpening
operations[46]) as ground-truths (GTs), and use them to synthesize the LR
images to build the LR-HR training pairs. In the early stage, bicubic
downsampling is commonly used to synthesize the LR images from their HR
counterparts [8, 9, 23, 33, 42, 56]. However, the ISR models trained on such
HR-LR pairs can hardly generalize to real-world images whose degradation
process is much more complex. Therefore, some researchers proposed to collect
HR-LR image pairs by using long-short camera focal lengths [4, 3]. While such
a degradation process is more reasonable than bicubic downsampling, it only
covers a small subspace of possible image degradations. Recently, researchers
[12, 20, 30, 51, 46, 50, 34, 32, 59] have proposed to shuffle or combine
different degradation factors, such as Gaussian/Poisson noise, (an-)isotropic
blur kernel, downsampling/upsampling, JPEG compression and so on, to
synthesize LR-HR image pairs, largely improving the generalization capability
of ISR models to real-world images.

Though great progress has been achieved, existing LR-HR training pair
generation schemes have several limitations. First, the original HR images are
used as the GTs to supervise the ISR model training. However, the perceptual
quality of HR images may not be high enough (Fig. 1 shows an example),
limiting the performance of the trained ISR models. Second, existing schemes
do not consider much human perception in GT generation, and the trained ISR
models tend to produce over-smoothed results. When the adversarial losses[27,
48, 40] are used to improve the ISR details, many unpleasant artifacts can be
introduced.

In order to tackle the aforementioned challenges, we propose a human guided GT
data generation strategy to train perceptually more realistic ISR (Real-ISR)
models. First, we elaborately train multiple image enhancement models to
improve the perceptual quality of HR images. Meanwhile, one LR image can have
multiple enhanced HR counterparts instead of only one. Second, to discriminate
the visual quality between the original and enhanced images, human subjects
are introduced to annotate the regions in enhanced HR images as “Positive”,
“Similar” or “Negative” samples, which represent better, similar or worse
perceptual quality compared with the original HR image. Consequently, a human
guided multiple-GT image dataset is constructed, which has both positive and
negative samples.

With the help of human annotation information in our dataset, positive and
negative LR-GT training pairs can be generated (examples of the positive and
negative GTs can be seen in Fig. 1), and a new loss function is proposed to
train the Real-ISR models. Extensive experiments are conducted to validate the
effectiveness and advantages of the proposed GT image generation strategy.
With the same backbone, the Real-ISR models trained on our dataset can produce
more perceptually realistic details with less artifacts than models trained on
the current datasets.

##  2 Related Work

According to how the LR-HR image pairs are created, the existing ISR methods
can be categorized into three major groups: simple degradation based, long-
short focal length based, and complex degradation based methods.

Simple Degradation based Training Pairs. Starting from SRCNN[8, 9], most of
the deep learning based ISR methods synthesize the LR images from their HR
counterparts by bicubic downsampling or direct downsampling after Gaussian
smoothing. By using such a simple degradation model to generate a large amount
of training data, researchers focus more on the ISR network module design,
such as residual[23]/dense[58] connection, channel-attention[56, 6, 17],
multiple receptive field[28, 16] or self-attention[5, 29, 54]. The fidelity
based measures, such as PSNR and SSIM[49], are used to evaluate and compare
the performance of different ISR methods. Later on, many works[27, 40, 47, 48,
41, 39, 35, 31] have been developed to adopt the Generative Adversarial
Network (GAN)[11] techniques to train Real-ISR models so as to produce photo-
realistic textures and details.

Long-short Focal Length based Training Pairs. Instead of synthesizing LR-HR
pairs using simple degradation operators, researchers have also tried to use
long-short camera focal length to collect real-world LR-HR pairs. The
representative works include CameraSR[4] and RealSR[3]. The former builds a
dataset using DSLR and mobile phone cameras to model degradation between the
image resolution and field-of-view. The latter utilizes different focal
lengths of the DSLR camera to shot the same scene at different resolutions,
and employs an image registration method to crop and align the LR-HR image
pairs. Nonetheless, ISR models trained on those datasets might fail when
applied to images from different resources (_e.g_., different degradation,
different focal length and cameras).

Figure 2: Illustration of our human guided ground-truth (GT) generation
process. We first train four image enhancement models to enhance the original
high-resolution (HR) image, and then extract the patches which have rich
textural and structural details while having certain differences between the
original and enhanced versions. Finally, human subjects are involved to
annotate the extract patches as “Positive”, “Similar” and “Negative” samples.

Complex Degradation based Training Pairs. The image degradation in real-world
scenarios can be too complex to model using a simple operator. To enable the
Real-ISR models having higher generalization capability, BSRGAN[51] and Real-
ESRGAN[46] have been proposed to synthesize LR-HR training pairs with more
complex image degradations. They employ a set of degradation factors, such as
different types of noise, blur kernels, scaling factors, JPEG compression,
_etc_., to enlarge the degradation space. BSRGAN[51] shuffles and combines
different degradations, while Real-ESRGAN[46] employs a two-stage synthesis
progress. In DASR[32], Liang _et al_. partitioned the complex degradation
space into different levels, and proposed a degradation adaptive method for
Real-ISR.

Other Training Pairs. Beside the above three groups of ISR methods,
MCinCGAN[57] and Pseudo-SR[36] utilize unpaired training images to do
unsupervised learning. They utilize one or more discriminators to tell the HR
GT from the unpaired SR output. AdaTarget[21] employs a transformation CNN
block to generate a training-friendly GT from the original GT during the
training progress. Nevertheless, the quality of the generated training-
friendly GT might not have a good perception quality.

##  3 Human Guided Ground-truth Generation

###  3.1 Overview

As discussed in Section 2, almost all existing methods [8, 9, 48, 15, 37, 55,
25, 34] directly take the HR images as the GT to construct the training pairs.
Unfortunately, the perceptual quality of many HR images may not be good enough
to serve as GTs, limiting the performance trained Real-ISR models. Therefore,
we propose to enhance the quality of HR images so that they can better serve
as GTs. In particular, human guidance can be introduced in the GT generation
process so that perceptually more realistic Real-ISR models can be trained.

As illustrated in Fig. 2, the proposed human guided GT generation method has
three steps. First, we elaborately train multiple image enhancement models to
improve the perceptual quality of HR images. Second, those patches which have
enough textural and structural details and have certain differences between
the enhanced version and the original version are extracted. Third, human
subjects are introduced to discriminate the visual quality between the
enhanced patches and the original patch, and label them as “Positive” (_i.e_.,
better quality), “Similar” (_i.e_., similar quality) or “Negative” (_i.e_.,
worse quality) samples. In the following subsections, we describe these three
steps in detail.

###  3.2 Design of the Enhancement Models

In order to generate visually more pleasing GTs from the original HR image, we
train multiple image enhancement models and apply them to the HR image. To
this end, the commonly used DF2K-OST dataset (including 800 high quality
images from DIV2K[1], 2650 high quality images from Flickr2K[43] and 10,324
images from OST[47]) is employed. The original images are denoted by $I^{H}$,
and the low quality ones, denoted by $I^{L}$, are degraded from $I^{H}$ by
using the following degradation model [51, 46]:

$I^{L}=[(I^{H}\otimes K)_{R}+V]_{J},$ (1)

where $K$ means isotropic/an-isotropic blur kernel, $R$ means resize
operation, $V$ is Gaussian/Poisson noise and $J$ denotes JPEG compression.
With ($I^{L}$, $I^{H}$) as training pairs, we can train enhancement models.
Note that before inputting the low-quality image $I^{L}$ into the model, we
resize it to the size of $I^{H}$ since here we are training enhancement
models, where the input and output have the same size.

Considering that the quality of HR image to be further enhanced is generally
not bad, we deliberately control the degradation settings in Eq. (1) within
weak to middle levels. Otherwise, the learned models can over-enhance the HR
images and generate many artifacts. Since the major issues of real world
images are noise corruption and blurring, we employ two degradation settings,
one focusing on processing slightly high noise and the other focusing on
dealing with slightly strong blur. The detailed degradation settings can be
found in the supplementary file.

We select one CNN-based network RCAN[56] and one transformer-based network
ELAN[54] as the backbones of our enhancer. RCAN[56] adopts deep residual
learning together with channel-attention[18], while ELAN[54] employs a multi-
scale self-attention[44] block to extract long-range independence. We remove
the up-sampling layer in those models since the input and output share the
same size in our case. We choose both CNN and transformer as backbones because
though transformers have stronger capability in restoring large scale
structures and repetitive patterns, CNN can better characterize some small
scale and local image details. With the two different degradation settings and
two different backbones, we train four image enhancement models with $L_{1}$,
perceptual and adversarial losses. The UNet discriminator[46] is used in
adversarial training.

###  3.3 Patch Selection and Annotation

We apply the trained four enhancement models to 1,600 HR images collected from
three representative resources: 1) 800 images from the DIV2K[1] dataset; 2)
400 images from Internet which could be used for free, such as Pixabay
(https://pixabay.com) and Unsplash (https://unsplash.com); 3) 400 images shot
by us using mobile phones. Note that though those HR images have high
resolution (2K$\sim$4K), they could contain certain noise, blurred details or
other real-world degradations, as we shown in Fig. 1. It is expected that
their perceptual quality can be improved by our enhancement models so that
they can better serve as GTs in Real-ISR model training.

After applying the four enhancement models to the 1,600 HR images, we obtain
6,400 enhanced images. However, it is inappropriate to directly take them as
GTs. On one hand, many regions in these images are smooth and less
informative. On the other hand, there is no guarantee that the enhancement
models can always produce perceptually better outputs in all regions.
Therefore, we extract patches from those enhanced images and invite human
volunteers to label them. In specific, we randomly crop $512*512$ patches from
each image with the overlapping area less than $1/2$ of patch area. We then
filter out the patches that have large smooth background regions according to
the quantity of details and textures, which is measured by the standard
deviation (std) of the patch in image domain and the std of high-frequency
components in a Laplacian pyramid. At last, we remove the patches on which the
difference between the original version and the enhanced version is small
(_i.e_., no much enhancement). The patch selection process avoids the cost of
annotating flat patches, and can speed up the training process since flat
patches have small gradients. Finally, we select 20,193 groups of patches of
$512*512$ size, each group having one original HR patch and 4 enhanced
patches.

We then invite 60 volunteers with different background to annotate the quality
of enhanced patches by comparing them with the original HR patch. A software
program, whose interface is shown in the supplementary file, is developed for
this purpose. The original patch is positioned at the left side of the screen,
while the four enhanced versions are located on the right side in random
order. Those patches whose perceptual quality is better than the original one
are labeled as “Positive”, and the patches with worse perceptual quality are
labeled as “Negative”. In case the quality is tied, the enhanced patch will be
labeled as “Similar”. Before annotating, all volunteers are briefly trained to
ensure that they will focus on the image perceptual quality (_e.g_.,
sharpness, noise, details, artifacts, _etc_.) but not on the image content.

###  3.4 Statistics of the Annotated Dataset

We invite 60 volunteers to annotate the 20,193 patch groups, each consisting
of an original HR patch and 4 enhanced patches. Each group is annotated by 3
different volunteers, and each volunteer is assigned with about 1,010 groups
to annotate. In total, we obtain 20,193 groups of annotations, and
$20,193\times 4\times 3=242,316$ annotated patches. The average annotation
time for one group is 22.79s.

Distribution of the patch annotations. Tab. 1 shows the distribution of
“Positive”, “Similar” and “Negative” labels for each enhancement model, as
well as the overall distribution. We see that there are overall 176,042
“Positive” (72.65%), 50,904 “Similar” (21.00%) and 15,370 “Negative” (6.35%)
patches. Such statistics imply that our enhancement models improve the visual
quality of HR patches in most cases, but there are indeed some bad cases.

Distribution of the final patch labels. For each of the $20,193\times
4=80,772$ enhanced patches, we have three annotations from three different
volunteers. We take the majority as the final label of the patch, _i.e_., if
one patch has two or three same annotations, it will be labeled by that
annotation. In case the three annotations are different from each other
(_i.e_., one “Positive”, one “Similar” and one “Negative”), the final label is
marked as “Similar”. Tab. 2 shows the distribution of the final labels of the
enhanced patches. We can see that finally there are 63,583 “Positive”
(78.72%), 14,675 “Similar” (18.17%) and 2,514 “Negative” (3.11%) patches. Most
of the final labels are “Positive” ones, and only a small portion (3.11%) are
“Negative” ones. The maximum divergence of “Positive” labels is 3,329 (5.24%)
between Model 2 and Model 3. The examples of “Positive”, “Similar” and
“Negative” patches can be found in the supplementary file.

Table 1: The distribution of annotations in our dataset. There are 20,193
groups of patches, while each group consists of an original HR patch and 4
enhanced patches. Each enhanced patch is annotated by 3 different volunteers,
resulting in a total of $20,193\times 4\times 3=242,316$ annotations.

Label | Enhance Model | Total
---|---|---
1 | 2 | 3 | 4
Positive | 42362 | 39031 | 47251 | 47398 | 176042
Similar | 14623 | 17615 | 10259 | 8407 | 50904
Negative | 3594 | 3933 | 3069 | 4774 | 15370
Total | 60579 | 60579 | 60579 | 60579 | 242316

Table 2: The distribution of final patch labels in our dataset. There are
$20,193\times 4=80,772$ enhanced patches, each having three annotations. We
take the majority annotation label as the final label of each patch.

Final Label | Enhance Model | Total
---|---|---
1 | 2 | 3 | 4
Positive | 15250 | 13907 | 17236 | 17190 | 63583
Similar | 4379 | 5635 | 2517 | 2144 | 14675
Negative | 564 | 651 | 440 | 859 | 2514
Total | 20193 | 20193 | 20193 | 20193 | 80772

Table 3: The distribution of the number ($0\sim 4$) of final “Positive”
patches per group in our dataset.

“Positive” Count | 0 | 1 | 2 | 3 | 4 | Total
---|---|---|---|---|---|---
Groups count | 1267 | 996 | 2616 | 3901 | 11413 | 20193

Distribution of the number of final “Positive” patches per group. For each
group of patches, there can be $0\sim 4$ final “Positive” samples. Tab. 3
shows the distribution of the number of final “Positive” patches per group.
One can see that among the 20,193 groups, 11,413 (56.52%) groups have 4
“Positive” patches, 3,901 (19.32%) have 3 “Positive” patches, 2,616 (12.95%)
have 2 “Positive” patches, 996 (4.93%) have 1 “Positive” patch, and 1,267
(6.28%) have none. We will use those “Positive” patches as “Positive” GTs, and
those “Negative” patches as “Negative” GTs to train Real-ISR models. The
patches with “Similar” labels are not employed in our training progress.

##  4 Training Strategies

As described in Sec. 3, for an original HR patch, denoted by $I^{H}$, we may
have several (less than 4) positive GTs, denoted by $I^{Pos}$, and several
negative GTs, denoted by $I^{Neg}$. To construct the positive or negative LR-
GT pairs for Real-ISR model training, we apply the degradation model in Eq. 1
to $I^{H}$ and obtain the corresponding LR image, denoted by $I^{L}$. (The
setting of degradation parameters will be discussed in Sec. 5.1). In total,
there are 63,583 positive LR-GT pairs $(I^{L},I^{Pos})$ and 2,514 negative LR-
GT pairs $(I^{L},I^{Neg})$. Note that in our dataset, one LR image may
correspond to multiple positive GTs or negative GTs.

Training with positive pairs only. By removing those groups that do not have
any positive GT from the 20,193 training groups, we have 18,926 groups with
$1\sim 4$ GTs, and 63,583 positive LR-GT pairs to train Real-ISR models. As in
previous works [51, 46], we employ the $L_{1}$ loss, perceptual loss $L_{p}$
and GAN loss $L_{GAN}$ to train the model. Since one LR image $I^{L}$ may have
multiple positive GTs, each time we randomly choose one positive GT to
calculate the $L_{1}$, $L_{p}$ and $L_{GAN}$ losses of the corresponding LR
image $I^{L}$, and update the discriminator and generator networks. The
overall training loss is as follows:

$L_{Total}=\alpha L_{1}+\beta L_{p}+\gamma L_{adv},$ (2)

where $\alpha$, $\beta$ and $\gamma$ are balance parameters.

Training with both positive and negative pairs. By filtering out those groups
that only contain “Similar” GTs, we obtain 19,272 groups that have at least
one “Positive” or “Negative” GT, totally 63,583 positive LR-GT pairs and 2,514
negative LR-GT pairs. When training with the positive GTs, we adopt the same
strategy as described above. For each negative LR-GT pair, we introduce a
negative loss, denoted by $L_{neg}$, to update the model.

It is observed that most of the negative GTs have over-sharpened details,
strong noise or false details (example images are provided in the
supplementary file). Inspired by LDL[31], we build a map
$\boldsymbol{M}^{Neg}$ to indicate the local residual variation of a negative
GT, which is defined as
$\boldsymbol{M}^{Neg}_{i,j}=var(\boldsymbol{R}^{Neg}_{i,j}(3,3))^{a}$, where
$\boldsymbol{R}^{Neg}=|I^{Neg}-I^{H}|$ is the residual between the original HR
image and the negative GT, $\boldsymbol{R}^{Neg}_{i,j}(3,3)$ is a local
$3\times 3$ window of $\boldsymbol{R}^{Neg}$ centered at $(i,j)$, $var$
denotes the variance operation and $a$ is the scaling factor (we set $a$ to
$\frac{3}{4}$ in our experiments).

Similarly, we can build a residual variation map
$\boldsymbol{M}^{Pos}_{i,j}=var(\boldsymbol{R}^{Pos}_{i,j}(3,3))^{a}$ for the
positive GT, where $\boldsymbol{R}^{Pos}=|I^{Pos}-I^{H}|$. At location
$(i,j)$, if the negative residual variation is higher than the positive one,
we identify this pixel at $I^{Neg}$ as a truly negative pixel, which should be
used to update the model. Therefore, we first define an indication map
$\boldsymbol{M}^{Ind}_{i,j}$:

$\boldsymbol{M}^{Ind}_{i,j}=\left\\{\begin{array}[]{rcl}0,&&{\boldsymbol{M}^{Neg}_{i,j}<=\boldsymbol{M}^{Pos}_{i,j}}\\\
\boldsymbol{M}^{Neg}_{i,j},&&{\boldsymbol{M}^{Neg}_{i,j}>\boldsymbol{M}^{Pos}_{i,j}}\end{array}\right.$
(3)

and then define the negative loss $L_{neg}$ as follows:

$L_{neg}=||\boldsymbol{M}^{Ind}\odot(I^{Neg}-I^{SR})||_{1},$ (4)

where $\odot$ means dot product.

Finally, the overall training loss is defined as:

$L_{Total}=\alpha L_{1}+\beta L_{p}+\gamma L_{adv}-\delta L_{neg},$ (5)

where $\alpha$, $\beta$, $\gamma$ and $\delta$ are balance parameters.

##  5 Experimental Results

###  5.1 Experiment Setup

To validate the effectiveness of our human guided GT (HGGT) dataset and the
role of negative GTs, we perform two sets of experiments. First, in Sec. 5.2,
we train several representative Real-ISR models, such as Real-ESRGAN[46],
BSRGAN[51], AdaTarget[21] and LDL[31] on the DF2K-OST dataset and our HGGT
dataset, and compare their performance. Second, in Sec. 5.3, we train two
commonly used Real-ISR backbones (RRDB[48, 46] and SwinIR[29]) on our dataset
by using only the postive GTs and using both the positive and negative GTs.

Implementation details. Before training a model on our dataset, we first pre-
train it on the DF2K-OST dataset by using the pixel-wise $\ell_{1}$ loss to
get a stable initialization. Since the original degradation settings in Real-
ESRGAN[46] and BSRGAN[51] is too strong to use in practical ISR applications,
we adopt a single-stage degradation process, including blur, noise, down-
sampling and JPEG compression with moderate intensities. Detailed settings and
visual examples are provided in the supplementary file. For the two backbones,
RRDB and SwinIR, we utilize the UNet discriminator[46] for adversarial
training, resulting in a RRDB-GAN model and a SwinIR-GAN model.

We conduct Real-ISR experiments with scaling factor 4 in this paper. We
randomly crop training patches of size $256*256$ from the GT images, and
resize the corresponding regions in the LR images to $64*64$. The batch size
is set to 12 for RRDB backbone and 8 for SwinIR backbone to save GPU memory.
We train our model on one NVIDIA RTX 3090 GPU for 300K iterations using the
Adam[24] optimizer. The initial learning rate is set to $1e-4$, and we halve
it after 200K iterations for RRDB backbone, and 200K, 250K, 275K and 287.5K
iterations for SwinIR backbone. The balance parameters $\alpha$, $\beta$,
$\gamma$ and $\delta$ in Eq. 5 are set to 1, 1, 0.1 and 300, respectively.
$\delta$ is set much larger than others because the number of negative GTs is
much smaller than positive ones.

Testing set. To evaluate the performance of Real-ISR models trained on our
dataset quantitatively, we construct a test set using the same steps as in the
construction of our training set. In specific, $100$ patch groups with at
least 2 ‘Positive’ GTs are constructed. The input LR patches are generated by
using the same degradation process as in the training process. The LR patches
together with their GTs are used to quantitatively evaluate the Real-ISR
models. We denote this dataset as Test-100.

Evaluation protocol. For the quantitative evaluation on Test-100, we adopt the
commonly used PSNR, SSIM [49] LPIPS [53] and DISTS[7] as quality metrics.
Since in Test-100 one LR image has at least 2 positive GTs, we average the
PSNR/SSIM/LPIPS/DISTS scores respectively over the multiple positive GTs as
the final scores. For the qualitative evaluation, we invite 12 volunteers to
perform subjective assessment, and report the count of preferred Real-ISR
models as the user study results.

Original&Positive GT

Real-ESRGAN

BSRGAN

LDL

AdaTarget-GAN

Figure 3: Visual comparison of state-of-the-art models trained on the DF2K-OST
and our proposed HGGT datasets. The 1st and 3rd rows show the results of
models trained on DF2K-OST, while the 2nd and 4th rows show the results of
models trained on ours positive GTs. The left column shows the original GT and
the positive GT in our dataset. Please zoom in for better observation.

Bicubic

Trained on Ori HR

Trained on Pos GT

Trained on Pos+Neg GT

Positive GT

Figure 4: Visualizations of RRDB-GAN ans SwinIR-GAN models trained on the
original HR (Ori HR) patches, positive GTs (Pos GT) only, and both positive
and negative GTs (Pos+Neg GT) in our HGGT dataset. The top and bottom rows
show the results of RRDB-GAN and SwinIR-GAN, respectively. From left to right
are the results of bicubic interpolation and the models trained on the Ori HR,
Pos GT, Pos+Neg GT, respectively. Please zoom in for better observation.

###  5.2 DF2K-OST Dataset vs. Our HGGT Dataset

We first evaluate the effectiveness of the proposed dataset by training
representative Real-ISR models respectively on the DF2K-OST dataset and the
positive GTs of our HGGT dataset. Four state-of-the-art Real-ISR models are
employed, including Real-ESRGAN[46], BSRGAN[51], AdaTarget[21] and LDL[31].
For Real-ESRGAN and BSRGAN, we adjust the degradation parameters so that the
quality of synthesized training LR images is comparable to the LR images in
our test set. For AdaTarget and LDL, we use the single-stage degradation as
explained in Sec. 5.1, and employ the loss functions in the original papers.
All models are firstly pre-trained on DF2K-OST with $\ell_{1}$ loss. The UNet
discriminator[46] is used for adversarial training in our experiments.
Quantitative comparison are reported in Table 4 and visual comparisons are
shown in Figure 3.

As shown in Table 4, training on our HGGT dataset leads to much better
LPIPS/DISTS scores against the DF2K-OST dataset. Specifically, the LPIPS/DISTS
scores are significantly improved by $10.14$%/$12.40$%, $16.30$%/$15.27$%,
$17.45$%/$18.91$% and $19.23$%/$21.99$%, respectively, for Real-ESRGAN,
BSRGAN, LDL and AdaTarget-GAN. This indicates a clear advantage of perceptual
quality brought by our dataset. Some visual examples are shown in Figure 3.
One can see that the models trained on our positive GTs can produce
perceptually more pleasing results against the models trained on DF2K-OST. The
reconstructed images by our dataset have sharper textures and richer details.
This is because the original GTs in the DF2K-OST dataset have mixed visual
qualities, where a large number of local patches are smooth. In comparison, in
our HGGT dataset, the perceptual quality of positive GTs is much enhanced, and
the smooth or artifactual patches are mannually removed. These improvements on
the training data bring clear advantages to the trained Real-ISR models. More
visual results are put in the supplementary file.

As a common problem of GAN-based models, the superior perceptual quality
sacrifices the pixel-wise fidelity which is depicted by PSNR and SSIM. This
trade-off, which is mainly caused by the ill-posed nature of the image
restoration tasks, has been discussed in previous researches[51]. It is well-
known that the pixel-wise metrics do not correlate well to the visual
quality[27, 40, 41]. In addition, in our proposed HGGT dataset, the perceptual
quality of GT is improved by using GAN-based enhancement models so that the
pixel-wise correlations may not be well-preserved in the data. However, human
observers generally prefer the enhanced images in our annotation process,
while the perceptually more pleasing results demonstrate the significance of
our proposed HGGT dataset in improving the upper bound of the Real-ISR tasks.

The main goal of the proposed HGGT dataset is to improve the perceptual
quality of Real-ISR outputs by introducing human perceptions into the training
pair generation. We perform a user study to validate the effectiveness of our
strategy by inviting 12 volunteers to evaluate the Real-ISR results on the
Test-100 dataset. For each of the four Real-ISR methods, _i.e_., Real-ESRGAN,
BSRGAN, AdaTarget-GAN and LDL, the two models trained on the DF2K-OST dataset
and the positive GTs of our HGGT dataset are compared. Each time, the Real-ISR
results of the two models on the same LR input are shown to the volunteers in
random order, and the volunteers are asked to chose the perceptually better
one based on their evaluation. The statistics of the user study are shown in
Fig. 5. It should be noted the volunteers invited in this user study do not
participate in the annotation process of our dataset.

As shown in Fig. 5, the majority of participants (more than 80% for all tests)
prefer the models trained on our HGGT dataset. This validates the
effectiveness of the proposed approach and the dataset, which can be plug-and-
play to most of the existing Real-ISR methods and improve their performance by
a large margin. For the images where models trained on DF2K-OST are selected,
we observe that they mostly contain much flat and smooth regions, and the
results of the two models are actually very close.

###  5.3 The Effectiveness of Negative GTs

We then evaluate the effectiveness of the negative GTs in our HGGT dataset. We
first train the baseline model on the original HR images that are used to
build our dataset. Then, we train the models on positive GTs only by using Eq.
(2), as illustrated in Section 4. Finally, we train the models on both
positive and negative GTs by using Eq. (5). The CNN-based RRDB and
transformer-based SwinIR backbones are used to train Real-ISR models. Due to
the limit of space, quantitative comparisons of the trained models are
reported in the supplementary file.

Table 4: The quantitative results of different Real-ISR models trained on
DF2K-OST and our HGGT datasets on Test-100.

Method |  | Train
---
Dataset
PSNR/SSIM/LPIPS/DISTS
Real-ESRGAN | DF2K-OST | 21.9797/0.6173/0.2593/0.1806
Positive GT | 21.5379/0.6078/0.2330/0.1582
BSRGAN | DF2K-OST | 21.7083/0.6092/0.2865/0.1880
Positive GT | 20.9037/0.5898/0.2398/0.1593
LDL | DF2K-OST | 22.4724/0.6394/0.2304/0.1676
Positive GT | 22.0190/0.6325/0.1902/0.1359
AdaTarget-GAN | DF2K-OST | 22.3944/0.6360/0.2335/0.1687
Positive GT | 21.9216/0.6301/0.1886/0.1316

Visual comparisons are shown in Figure 4, which provides more intuitive
evidences on the effectiveness of the annotated negative GTs. As shown in the
second column, the models trained on original HR images yield blurry details
and irregular patterns, especially on the area with dense textures. This is
mainly caused by the low and mixed quality of the original HR images. In
contrast, training on our positive GTs can produce much sharper textures and
richer details, whereas there remain some false details and visual artifacts
(see the windows of the building). Further, training on both positive and
negative GTs leads to a more balanced visual performance. Some over-enhanced
local pixels can be suppressed, while the textures remain sharp and regular.
This is owing to the effective annotation of negative samples, which bring
useful human perception guidance into the data for model training. More visual
results can be found in the supplementary file.

Figure 5: User study results on the Real-ISR models trained on the DF2K-OST
dataset (the blue bar) and the positive GTs in our HGGT dataset (the red bar).

##  6 Conclusion

In this paper, we elaborately designed a human guided ground-truth (GT)
generation method for realistic image super-resolution (Real-ISR). We first
trained four image enhancement models to improve the perceptual quality of
original high resolution images, and then extracted structural and textural
patches from the enhanced images. Finally, human subjects were invited to
annotate the perceptual quality of extracted patches as positive and negative
GTs, resulting in the human guided ground-truth (HGGT) dataset. The sharper
textures and richer details in the positive GTs could largely improve the
performance of trained Real-ISR models, while the negative GTs could provide
further guidance for the model to avoid generating visual artifacts. Extensive
experiments validated the effectiveness of the proposed HGGT dataset and the
training strategies.

Acknowledgement. We thank Dr. Lida LI for providing support in GPU server
configuration and the many people participating in data collection and
annotation.

## References

  * [1] Eirikur Agustsson and Radu Timofte.  Ntire 2017 challenge on single image super-resolution: Dataset and study.  In CVPRW, pages 126–135, 2017.
  * [2] Sefi Bell-Kligler, Assaf Shocher, and Michal Irani.  Blind super-resolution kernel estimation using an internal-gan.  NeurIPS, 32, 2019.
  * [3] Jianrui Cai, Hui Zeng, Hongwei Yong, Zisheng Cao, and Lei Zhang.  Toward real-world single image super-resolution: A new benchmark and a new model.  In ICCV, pages 3086–3095, 2019.
  * [4] Chang Chen, Zhiwei Xiong, Xinmei Tian, Zheng-Jun Zha, and Feng Wu.  Camera lens super-resolution.  In CVPR, pages 1652–1660, 2019.
  * [5] Hanting Chen, Yunhe Wang, Tianyu Guo, Chang Xu, Yiping Deng, Zhenhua Liu, Siwei Ma, Chunjing Xu, Chao Xu, and Wen Gao.  Pre-trained image processing transformer.  In CVPR, pages 12299–12310, 2021.
  * [6] Tao Dai, Jianrui Cai, Yongbing Zhang, Shu-Tao Xia, and Lei Zhang.  Second-order attention network for single image super-resolution.  In CVPR, pages 11065–11074, 2019.
  * [7] Keyan Ding, Kede Ma, Shiqi Wang, and Eero P Simoncelli.  Image quality assessment: Unifying structure and texture similarity.  IEEE TPAMI, 44(5):2567–2581, 2020.
  * [8] Chao Dong, Chen Change Loy, Kaiming He, and Xiaoou Tang.  Learning a deep convolutional network for image super-resolution.  In ECCV, pages 184–199. Springer, 2014.
  * [9] Chao Dong, Chen Change Loy, Kaiming He, and Xiaoou Tang.  Image super-resolution using deep convolutional networks.  IEEE TPAMI, 38(2):295–307, 2015.
  * [10] Chao Dong, Chen Change Loy, and Xiaoou Tang.  Accelerating the super-resolution convolutional neural network.  In ECCV, pages 391–407. Springer, 2016.
  * [11] Ian Goodfellow, Jean Pouget-Abadie, Mehdi Mirza, Bing Xu, David Warde-Farley, Sherjil Ozair, Aaron Courville, and Yoshua Bengio.  Generative adversarial networks.  Communications of the ACM, 63(11):139–144, 2020.
  * [12] Jinjin Gu, Hannan Lu, Wangmeng Zuo, and Chao Dong.  Blind super-resolution with iterative kernel correction.  In CVPR, pages 1604–1613, 2019.
  * [13] Muhammad Haris, Greg Shakhnarovich, and Norimichi Ukita.  Deep back-projectinetworks for single image super-resolution.  IEEE TPAMI, 43(12):4323–4337, 2020.
  * [14] Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun.  Deep residual learning for image recognition.  In CVPR, pages 770–778, 2016.
  * [15] Xiangyu He, Zitao Mo, Peisong Wang, Yang Liu, Mingyuan Yang, and Jian Cheng.  Ode-inspired network design for single image super-resolution.  In CVPR, pages 1732–1741, 2019.
  * [16] Zewei He, Yanpeng Cao, Lei Du, Baobei Xu, Jiangxin Yang, Yanlong Cao, Siliang Tang, and Yueting Zhuang.  Mrfn: Multi-receptive-field network for fast and accurate single image super-resolution.  IEEE TMM, 22(4):1042–1054, 2019.
  * [17] Zewei He, Du Chen, Yanpeng Cao, Jiangxin Yang, Yanlong Cao, Xin Li, Siliang Tang, Yueting Zhuang, and Zheming Lu.  Single image super-resolution based on progressive fusion of orientation-aware features.  PR, page 109038, 2022.
  * [18] Jie Hu, Li Shen, and Gang Sun.  Squeeze-and-excitation networks.  In CVPR, pages 7132–7141, 2018.
  * [19] Gao Huang, Zhuang Liu, Laurens Van Der Maaten, and Kilian Q Weinberger.  Densely connected convolutional networks.  In CVPR, pages 4700–4708, 2017.
  * [20] Yan Huang, Shang Li, Liang Wang, Tieniu Tan, et al.  Unfolding the alternating optimization for blind super resolution.  NeurIPS, 33:5632–5643, 2020.
  * [21] Younghyun Jo, Seoung Wug Oh, Peter Vajda, and Seon Joo Kim.  Tackling the ill-posedness of super-resolution through adaptive target generation.  In CVPR, pages 16236–16245, 2021.
  * [22] Justin Johnson, Alexandre Alahi, and Li Fei-Fei.  Perceptual losses for real-time style transfer and super-resolution.  In ECCV, pages 694–711. Springer, 2016.
  * [23] Jiwon Kim, Jung Kwon Lee, and Kyoung Mu Lee.  Accurate image super-resolution using very deep convolutional networks.  In CVPR, pages 1646–1654, 2016.
  * [24] Diederik P Kingma and Jimmy Ba.  Adam: A method for stochastic optimization.  arXiv preprint arXiv:1412.6980, 2014.
  * [25] Xiangtao Kong, Hengyuan Zhao, Yu Qiao, and Chao Dong.  Classsr: A general framework to accelerate super-resolution networks by data characteristic.  In CVPR, pages 12016–12025, 2021.
  * [26] Wei-Sheng Lai, Jia-Bin Huang, Narendra Ahuja, and Ming-Hsuan Yang.  Deep laplacian pyramid networks for fast and accurate super-resolution.  In CVPR, pages 624–632, 2017.
  * [27] Christian Ledig, Lucas Theis, Ferenc Huszár, Jose Caballero, Andrew Cunningham, Alejandro Acosta, Andrew Aitken, Alykhan Tejani, Johannes Totz, Zehan Wang, et al.  Photo-realistic single image super-resolution using a generative adversarial network.  In CVPR, pages 4681–4690, 2017.
  * [28] Juncheng Li, Faming Fang, Kangfu Mei, and Guixu Zhang.  Multi-scale residual network for image super-resolution.  In ECCV, pages 517–532, 2018.
  * [29] Jingyun Liang, Jiezhang Cao, Guolei Sun, Kai Zhang, Luc Van Gool, and Radu Timofte.  Swinir: Image restoration using swin transformer.  In ICCVW, pages 1833–1844, 2021.
  * [30] Jingyun Liang, Guolei Sun, Kai Zhang, Luc Van Gool, and Radu Timofte.  Mutual affine network for spatially variant kernel estimation in blind image super-resolution.  In ICCV, pages 4096–4105, 2021.
  * [31] Jie Liang, Hui Zeng, and Lei Zhang.  Details or artifacts: A locally discriminative learning approach to realistic image super-resolution.  In CVPR, pages 5657–5666, 2022.
  * [32] Jie Liang, Hui Zeng, and Lei Zhang.  Efficient and degradation-adaptive network for real-world image super-resolution.  In ECCV, pages 574–591, 2022.
  * [33] Bee Lim, Sanghyun Son, Heewon Kim, Seungjun Nah, and Kyoung Mu Lee.  Enhanced deep residual networks for single image super-resolution.  In CVPRW, pages 136–144, 2017.
  * [34] Zhengxiong Luo, Yan Huang, Shang Li, Liang Wang, and Tieniu Tan.  Learning the degradation distribution for blind image super-resolution.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 6063–6072, 2022.
  * [35] Cheng Ma, Yongming Rao, Yean Cheng, Ce Chen, Jiwen Lu, and Jie Zhou.  Structure-preserving super resolution with gradient guidance.  In CVPR, pages 7769–7778, 2020.
  * [36] Shunta Maeda.  Unpaired image super-resolution using pseudo-supervision.  In CVPR, pages 291–300, 2020.
  * [37] Salma Abdel Magid, Yulun Zhang, Donglai Wei, Won-Dong Jang, Zudi Lin, Yun Fu, and Hanspeter Pfister.  Dynamic high-pass filtering and multi-spectral attention for image super-resolution.  In ICCV, pages 4288–4297, 2021.
  * [38] Yiqun Mei, Yuchen Fan, and Yuqian Zhou.  Image super-resolution with non-local sparse attention.  In CVPR, pages 3517–3526, 2021.
  * [39] Mohammad Saeed Rad, Behzad Bozorgtabar, Urs-Viktor Marti, Max Basler, Hazim Kemal Ekenel, and Jean-Philippe Thiran.  Srobb: Targeted perceptual loss for single image super-resolution.  In ICCV, pages 2710–2719, 2019.
  * [40] Mehdi SM Sajjadi, Bernhard Scholkopf, and Michael Hirsch.  Enhancenet: Single image super-resolution through automated texture synthesis.  In ICCV, pages 4491–4500, 2017.
  * [41] Jae Woong Soh, Gu Yong Park, Junho Jo, and Nam Ik Cho.  Natural and realistic single image super-resolution with explicit natural manifold discrimination.  In CVPR, pages 8122–8131, 2019.
  * [42] Ying Tai, Jian Yang, and Xiaoming Liu.  Image super-resolution via deep recursive residual network.  In CVPR, pages 3147–3155, 2017.
  * [43] Radu Timofte, Eirikur Agustsson, Luc Van Gool, Ming-Hsuan Yang, and Lei Zhang.  Ntire 2017 challenge on single image super-resolution: Methods and results.  In CVPRW, pages 114–125, 2017.
  * [44] Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N Gomez, Łukasz Kaiser, and Illia Polosukhin.  Attention is all you need.  NeurIPS, 30, 2017.
  * [45] Longguang Wang, Yingqian Wang, Xiaoyu Dong, Qingyu Xu, Jungang Yang, Wei An, and Yulan Guo.  Unsupervised degradation representation learning for blind super-resolution.  In CVPR, pages 10581–10590, 2021.
  * [46] Xintao Wang, Liangbin Xie, Chao Dong, and Ying Shan.  Real-esrgan: Training real-world blind super-resolution with pure synthetic data.  In ICCVW, pages 1905–1914, 2021.
  * [47] Xintao Wang, Ke Yu, Chao Dong, and Chen Change Loy.  Recovering realistic texture in image super-resolution by deep spatial feature transform.  In CVPR, pages 606–615, 2018.
  * [48] Xintao Wang, Ke Yu, Shixiang Wu, Jinjin Gu, Yihao Liu, Chao Dong, Yu Qiao, and Chen Change Loy.  Esrgan: Enhanced super-resolution generative adversarial networks.  In ECCVW, pages 0–0, 2018.
  * [49] Zhou Wang, Alan C Bovik, Hamid R Sheikh, and Eero P Simoncelli.  Image quality assessment: From error visibility to structural similarity.  IEEE TIP, 13(4):600–612, 2004.
  * [50] Zongsheng Yue, Qian Zhao, Jianwen Xie, Lei Zhang, Deyu Meng, and Kwan-Yee K Wong.  Blind image super-resolution with elaborate degradation modeling on noise and kernel.  In CVPR, pages 2128–2138, 2022.
  * [51] Kai Zhang, Jingyun Liang, Luc Van Gool, and Radu Timofte.  Designing a practical degradation model for deep blind image super-resolution.  In ICCV, pages 4791–4800, 2021.
  * [52] Kai Zhang, Wangmeng Zuo, and Lei Zhang.  Learning a single convolutional super-resolution network for multiple degradations.  In CVPR, pages 3262–3271, 2018.
  * [53] Richard Zhang, Phillip Isola, Alexei A Efros, Eli Shechtman, and Oliver Wang.  The unreasonable effectiveness of deep features as a perceptual metric.  In CVPR, pages 586–595, 2018.
  * [54] Xindong Zhang, Hui Zeng, Shi Guo, and Lei Zhang.  Efficient long-range attention network for image super-resolution.  In ECCV, pages 649–667, 2022.
  * [55] Xindong Zhang, Hui Zeng, and Lei Zhang.  Edge-oriented convolution block for real-time super resolution on mobile devices.  In ACM MM, pages 4034–4043, 2021.
  * [56] Yulun Zhang, Kunpeng Li, Kai Li, Lichen Wang, Bineng Zhong, and Yun Fu.  Image super-resolution using very deep residual channel attention networks.  In ECCV, pages 286–301, 2018.
  * [57] Yongbing Zhang, Siyuan Liu, Chao Dong, Xinfeng Zhang, and Yuan Yuan.  Multiple cycle-in-cycle generative adversarial networks for unsupervised image super-resolution.  IEEE TIP, 29:1101–1112, 2019.
  * [58] Yulun Zhang, Yapeng Tian, Yu Kong, Bineng Zhong, and Yun Fu.  Residual dense network for image super-resolution.  In CVPR, pages 2472–2481, 2018.
  * [59] Hongyi Zheng, Hongwei Yong, and Lei Zhang.  Unfolded deep kernel estimation for blind image super-resolution.  In ECCV, pages 502–518, 2022.


# RECURSIVE APPROACH TO THE MATRIX MOMENT PROBLEM

Abderrazzak Ech-charyfy<EMAIL_ADDRESS>Kaissar Idrissi
<EMAIL_ADDRESS>El Hassan Zerouali<EMAIL_ADDRESS>
###### Abstract

We investigate in this paper the truncated matrix moment problems of
Hamburger, Stieltjes and Hausdorff type by the mean of the description of non
negative matrix polynomials and by the use of linear recurrence relations. We
extend several known result in the scalar moment problem to the matrix moment
setting.

##  1 Introduction

The real and complex one-dimensional moment problems are well-known in
classical analysis and have been widely studied by many mathematicians and
engineers since the end of the $19th$ century. Stieltjes studied the problem
on the positive real line in [19], Hamburger on the real line in [10] and
Hausdorff on the bounded segment in [11]. The key point in the resolution of
the full moment problem is the use of Riesz’s theorem [17] on one hand and an
adequate description of positive polynomials on the other hand. The solution
of the one-dimensional moment problem is then given in terms of the positivity
of some associated Hankel matrices. The case of the truncated case is more
constructive and rely on some techniques from finite dimension linear algebra
together with extension approaches.

The the study of the full and the truncated moment problems are closely
related. In fact, J. Stochel in [20], showed that the solution of the full
problem is equivalent to the solution of the corresponding truncated problems.
For the truncated moment problems, R. Curto and L. Fialkow leaded a long and
complete study of the problem based a flatness concept, see for example
[4],[5] and [6], D. Kimsey and H. Woerdeman gave a solution based on
commutativity conditions of certain matrices [12]. There are also several
known mathematicians that tackled to problem using different approaches. K.
Schmüdgen uses the operator approach [18], F. H. Vasilescu the dimension
stability [21], C. Berg, and Y. Stochel, among others.

We adopt in this paper the recursive relation approach initiated by M. Rachidi
et all [1], [2] and we will extend the results from the paper [1] to the
matrix case.

This paper is divided into two parts. In the first part, we will recover in a
simple way some necessary and sufficient conditions for the one-dimensional
Stieltjes, Hamburger and Hausdorff matrix moment problems [15] to have a
solution. To do so, we will use the generalization of Riesz’s theorem in the
matrix setting that has been obtained by J. Cimprič and A. Zalar in [3], and
the characterization of matrix polynomials with real variables. In the second
part, we will extend the recurrence approach which was studied by R. BenTaher,
M. Rachidi and E. H. Zerouali [1] to the matrix case. This allows to solve the
corresponding truncated matrix moment problems.

##  2 Preliminaries

We denote by $\mathcal{M}_{p}(\mathbb{R})$ for the set of $p\times p$ matrices
with real entries and $\mathcal{M}_{p}(\mathbb{R}[X])$ the algebra of $p\times
p$ matrices with entries in $\mathbb{R}[X]$. Alternatively,
$\mathcal{M}_{p}(\mathbb{R}[X])$ is the algebra of polynomials with matrix
coefficients in $\mathcal{M}_{p}(\mathbb{R})$. That is every polynomial
$P\in\mathcal{M}_{p}(\mathbb{R}[X])$ has the form,

$P(X)=\displaystyle\sum_{k=0}^{n}A_{k}X^{k},\quad
A_{k}\in\mathcal{M}_{p}(\mathbb{R}).$

The following definition of the matrix-valued $K$-moment problem and and the
tracial $K$-moment problem can be found in [16]. Let
$S=(S_{n})_{n\in\mathbb{N}}$ be a sequence of symmetric matrices of order $p$
and $K\subset\mathbb{R}$ be a closed set. We consider the real valued function
$\mathbf{L}_{S}$ on $\mathcal{M}_{p}(\mathbb{R}[X])$ defined by:

$\mathbf{L}_{S}:\mathcal{M}_{p}(\mathbb{R}[X])\ni
P\equiv\displaystyle\sum_{k=0}^{m}A_{k}X^{k}\longmapsto\mathbf{L}_{S}(P):=\displaystyle\sum_{k=0}^{m}Tr\left(A_{k}S_{k}\right)\in\mathbb{R}.$

The linear form $\mathbf{L}_{S}$ is known as the Riesz functional associated
with the sequence $S$. The tracial $K$-moment problem is then introduced as
follows

###### Problem 1

The tracial $K$-moment problem asks when does there exist a positive semi-
definite matrix-valued measure $\sigma=(\sigma_{ij})_{1\leq i,j\leq p}$
satisfying the conditions:

$Supp(\sigma):=\displaystyle\bigcup_{1\leq i,j\leq p}Supp(\sigma_{ij})\subset
K,\quad\mathbf{L}_{S}(P)=\displaystyle\int_{K}Tr\big{(}P(x)d\sigma(x)\big{)}.$

In this case $\mathbf{L}_{S}$ is called a tracial $K$-moment functional and
$\sigma$ is its representing measure.

We also consider the matrix $K$-moment problem given next.

###### Problem 2

The matrix $K$-moment problem asks when does there exist a positive semi-
definite matrix-valued measure $\sigma=(\sigma_{ij})_{1\leq i,j\leq p}$
satisfying the conditions:

$Supp(\sigma):=\displaystyle\bigcup_{1\leq i,j\leq p}Supp(\sigma_{ij})\subset
K,\quad S_{k}=\displaystyle\int_{K}x^{k}d\sigma(x).$

We keep the usual definitions in the case of matrix moment problems. The
problem is called the Hamburger matrix moment problem if $K=\mathbb{R}$, the
Stieltjes matrix moment problem if $K=\left[0,+\infty\right[$ and the
Hausdorff matrix moment problem if $K=[0,1]$.

Recall the next two known results from [16] to be used in this paper.

######  Proposition 1 ([16, Proposition 2.8])

The following are equivalent:

  1. 1.

$S$ is a matrix-valued $K$-moment sequence;

  2. 2.

$\mathbf{L}_{S}$ is a tracial $K$-moment functional.

□

######  Theorem 1 ([16, Theorem 2.9])

The following are equivalent:

  1. 1.

$\mathbf{L}_{S}$ is a tracial $K$-moment functional;

  2. 2.

$\mathbf{L}_{S}(P)\geq 0$ for all $P\geq 0$ on $K$.

□

We will take advantage of these results to give new proofs for the solution of
the Hamburger, the Stieltjes and the Hausdorff matrix moment problems. Relying
on the characterization of positive matrix polynomials in each case, we are
able to simplify the known proofs in the two first cases and to refine the
result in the Hausdorff case.

##  3 Solutions of classical matrix moment problems.

###  3.1 Hamburger matrix moment sequences.

M.G. Krein introduced the Hamburger matrix moment problem in [15] in $1970$
where he also provided a solution relying on the extension of a symmetric
operator in an operator approach. We will give in this subsection a new proof
of Krein’s result. To this aim, we will use the next generalization of the
Fejér-Riesz theorem from [8], to recover the characterization of positive
matrix polynomials on $\mathbb{R}$.

######  Theorem 2 ([8, Theorem 2.1])

Let $L(z)=\displaystyle\sum_{k=-m}^{m}B_{k}z^{k}$ be a matrix Laurent
polynomial with coefficients in $\mathcal{M}_{p}(\mathbb{C})$, if $L(\xi)\geq
0$ for all $\xi\in\mathrm{T}=\left\\{z\in\mathbb{C}:|z|=1\right\\}$, then

$L(\xi)=R(\xi)^{*}R(\xi),\quad\xi\in\mathrm{T},$ (1)

for some matrix polynomial $R(z)=\displaystyle\sum_{k=0}^{m}R_{k}z^{k}$. □

We derive the next known lemma in a simple way.

###### Lemma 1

Let $P\in\mathcal{M}_{p}(\mathbb{R}[X])$. The following are equivalent:

  1. 1.

$P\geq 0$ on $\mathbb{R}$;

  2. 2.

$P$ is a sum of two symmetric squares of polynomials from
$\mathcal{M}_{p}\left(\mathbb{R}[X]\right)$.

□

###### Proof

Since $(ii)\Rightarrow(i)$ is trivial, it suffices to show
$(i)\Rightarrow(ii)$. If
$P(X)=\displaystyle\sum_{k=0}^{n}A_{k}X^{k}\in\mathcal{M}_{p}(\mathbb{R}[X])$
is such that $P(x)\geq 0$ for every $x\in\mathbb{R}$, the clearly $n=2m$ for
some integer $m$. Now the mapping
$\left]-\frac{\pi}{2},\frac{\pi}{2}\right[\ni t\longmapsto
tan(t)\in\mathbb{R}$ is bijective and writing $x=\tan(t)$ for some
$t\in\left]-\frac{\pi}{2},\frac{\pi}{2}\right[$ will give,

$P(x)=P(tan(t))=\displaystyle\frac{1}{\cos{(t)}^{2m}}\sum_{k=0}^{2m}A_{k}\sin{(t)}^{k}\cos{(t)}^{2m-k}.$

The function $Q(\sin{(t)},\cos{(t)})=\cos{(t)}^{2m}P\left(tan(t)\right)$
satisfies,

$Q(u,v)=v^{2m}P\left(\frac{u}{v}\right)\geq 0,\quad
u\in\mathbb{R},v\in\mathbb{R}^{*}.$

In particular

$Q(x,1)=P(x)\;\mbox{ for every }\;x\in\mathbb{R}.$

For $0\leq k\leq 2m$, we have:

$\sin{(x)}^{k}\cos{(x)}^{2m-k}=\sum_{p=0}^{k}\sum_{q=0}^{2m-k}\frac{1}{i^{k}2^{2m}}C_{p}^{k}C_{q}^{2m-k}e^{2ix(p+q-m)}.$

Moreover $0\leq p\leq k$ and $0\leq q\leq 2m-k$, hence $-m\leq p+q-m\leq m.$
Then there exist a matrix Laurent polynomial
$L(z)=\displaystyle\sum_{k=-m}^{m}B_{k}z^{k}$,
$B_{k}\in\mathcal{M}_{p}(\mathbb{C})$ such that:

$Q(\sin{(t)},\cos{(t)})=\sum_{k=0}^{2m}A_{k}\sin{(t)}^{k}\cos{(t)}^{2m-k}=L(e^{2it}).$
(2)

Since moreover, $L(\xi)\geq 0$ for all
$\xi\in\mathrm{T}=\\{z\in\mathbb{C}:|z|=1\\}$, according to Fejér-Riesz’s
theorem 2, we obtain

$L(e^{it})=R(e^{it})^{*}R(e^{it}),\quad t\in\mathbb{R},$ (3)

then

$Q(\sin{(t)},\cos{(t)})=R(e^{2it})^{*}R(e^{2it}).$

Now writing

$R(e^{2it})=A(\sin{(t)},\cos{(t)})+iB(\sin{(t)},\cos{(t)}).$ (4)

According to (2), (3) and (4), we get

$\displaystyle Q\left(\sin{(t)},\cos{(t)}\right)$ $\displaystyle=$
$\displaystyle{{}^{\top}}A(\sin{(t)},\cos{(t)})A(\sin{(t)},\cos{(t)})$
$\displaystyle+$ $\displaystyle
i{{}^{\top}}A(\sin{(t)},\cos{(t)})B(\sin{(t)},\cos{(t)})$ $\displaystyle-$
$\displaystyle i{{}^{\top}}B(\sin{(t)},\cos{(t)})A(\sin{(t)},\cos{(t)})$
$\displaystyle+$
$\displaystyle{{}^{\top}}B(\sin{(t)},\cos{(t)})B(\sin{(t)},\cos{(t))},$

where ${{}^{\top}M}$ stands for the usual transpose a matrix $M$. Since $Q$
has real coefficients, we deduce that,

${{}^{\top}}A(\sin{(t)},\cos{(t)}B(\sin{(t)},\cos{(t)})-{{}^{\top}}B(\sin{(t)},\cos{(t)})A(\sin{(t)},\cos{(t)})=0$

and then

$Q\left(\sin{(t)},\cos{(t)}\right)={{}^{\top}}A(\sin{(t)},\cos{(t)})A(\sin{(t)},\cos{(t)})+{{}^{\top}}B(\sin{(t)},\cos{(t)})B(\sin{(t)},\cos{(t)}).$

Finally

$P(x)={{}^{\top}}A(x,1)A(x,1)+{{}^{\top}}B(x,1)B(x,1).$

■

In sequel, denote $H_{m}=(S_{i+j})_{1\leq i,j\leq m},\mbox{
}EH_{m}=(S_{i+j+1})_{1\leq i,j\leq m}\mbox{ and }E^{2}H_{m}=(S_{i+j+2})_{1\leq
i,j\leq m}$. As linear forms $H_{m},EH_{m}\mbox{ and }E^{2}H_{m}$ are non
negative respectively if and only if $\displaystyle\sum_{0\leq i,j\leq
m}c_{i}S_{i+j}{{}^{\top}}c_{j}\geq 0$, $\displaystyle\sum_{0\leq i,j\leq
m}c_{i}S_{i+j+1}{{}^{\top}}c_{j}\geq 0$ and $\displaystyle\sum_{0\leq i,j\leq
m}c_{i}S_{i+j+2}{{}^{\top}}c_{j}\geq 0$ for all
$c_{0},c_{1},...,c_{m}\in\mathbb{R}^{p}$ respectively.

We give now a new proof of Krein theorem.

######  Theorem 3 ([15, Krein Theorem])

The following are equivalent:

  1. 1.

$S=(S_{n})_{n\in\mathbb{N}}$ is a Hamburger matrix moment sequence;

  2. 2.

$H_{m}$ is positive semi-definite for every $m\in\mathbb{N}$.

□

###### Proof

The assertion $(i)\Rightarrow(ii)$ is clear. For the converse, suppose that
$H_{m}$ is positive semi-definite for every $m\in\mathbb{N}$, and let us show
that $S$ is a Hamburger matrix moment sequence. From Proposition 1, this is
equivalent to show that $\mathbf{L}_{S}$ is a tracial moment functional on
$\mathbb{R}$ and using Theorem 1, it suffices to show that
$\mathbf{L}_{S}(P)\geq 0$ for all $P\geq 0$ on $\mathbb{R}$. Now according to
Lemma 1 it is enough to show that $\mathbf{L}_{S}({{}^{\top}}PP)\geq 0$ for
all $P\in\mathcal{M}_{p}(\mathbb{R}[X])$. To show the last assumption, let
$P(X)=\displaystyle\sum_{i=0}^{m}A_{i}X^{i}\in\mathcal{M}_{p}(\mathbb{R}[X])$,
we have
${{}^{\top}}PP(X)=\displaystyle\sum_{i,j=0}^{m}{{}^{\top}}A_{i}X^{i+j}A_{j}.$
Since
$Tr({{}^{\top}}A_{i}S_{i+j}A_{j})=\displaystyle\sum_{k=0}^{m}{{}^{\top}}(A_{j}e_{k})S_{i+j}A_{i}e_{k}\geq
0$, we obtain

$\mathbf{L}_{S}({{}^{\top}}PP)=\sum_{i,j=0}^{m}Tr({{}^{\top}}A_{i}S_{i+j}A_{j})=\displaystyle\sum_{k=0}^{m}\left(\sum_{i,j=0}^{m}{{}^{\top}}(A_{i}e_{k})S_{i+j}A_{j}e_{k}\right)\geq
0.$

■

###  3.2 Stieltjes matrix moment sequences.

The Stieltjes matrix moment problem was introduced by M.G. Krein and M.A.
Krasnosel’skiy [14]. They used the operator approach to obtain necessary and
sufficient conditions for the solvability of this problem.

We first show in a simple manner the next known lemma that describes positive
matrix polynomials on $[0,+\infty[$ needed in the proof of M.G. Krein and M.A.
Krasnosel’skiy result’s.

###### Lemma 2

Let $P\in\mathcal{M}_{p}(\mathbb{R}[X])$. The following are equivalent:

  1. 1.

$P\geq 0$ on $\left[0,+\infty\right[$;

  2. 2.

$P(X)=A(X)+XB(X)$, where $A$ and $B$ are sums symmetric of squares of
polynomials.

□

###### Proof

It suffices to show $(i)\Rightarrow(ii)$. For $P\geq 0$ on $[0,+\infty[$, the
polynomial $\displaystyle Q(y)=P(y^{2})$ in non negative on $\mathbb{R}$, and
from Lemma 1, we have

$\displaystyle
Q(y)=P(y^{2})={{}^{\top}}P_{1}(y)P_{1}(y)+{{}^{\top}}P_{2}(y)P_{2}(y).$

Now splitting odd and even degrees, we write

$Q(y)={{}^{\top}}(A_{1}(y^{2})+yA_{2}(y^{2}))(A_{1}(y^{2})+yA_{2}(y^{2}))+{{}^{\top}}(B_{1}(y^{2})+yB_{2}(y^{2}))(B_{1}(y^{2})+yB_{2}(y^{2})).$

It follows that

$\displaystyle Q(y)$
$\displaystyle={{}^{\top}}A_{1}(y^{2})A_{1}(y^{2})+y^{2}{{}^{\top}}A_{2}(y^{2})A_{2}(y^{2})+{{}^{\top}}B_{1}(y^{2})B_{1}(y^{2})$
$\displaystyle+y^{2}{{}^{\top}}B_{2}(y^{2})B_{2}(y^{2})+y({{}^{\top}}A_{1}(y^{2})A_{2}(y^{2})+{{}^{\top}}A_{2}(y^{2})A_{1}(y^{2})$
$\displaystyle+{{}^{\top}}B_{1}(y^{2})B_{2}(y^{2})+{{}^{\top}}B_{2}(y^{2})B_{1}(y^{2})).$

Since $Q(y)=Q(-y)$, we derive

$\displaystyle
Q(y)=\frac{1}{2}\left[{{}^{\top}}A_{1}(y^{2})A_{1}(y^{2})+{{}^{\top}}B_{1}(y^{2})B_{1}(y^{2})+y^{2}\big{(}{{}^{\top}}A_{2}(y^{2})A_{2}(y^{2})+{{}^{\top}}B_{2}(y^{2})B_{2}(y^{2})\big{)}\right].$

The change of variables $x=y^{2}$ gives

$P(x)=A(x)+xB(x),$

where
$A(x)=\displaystyle\frac{1}{2}\big{[}{{}^{\top}}A_{1}(x)A_{1}(x)+{{}^{\top}}B_{1}(x)B_{1}(x)\big{]}$
and
$B(x)=\displaystyle\frac{1}{2}\big{[}{{}^{\top}}A_{2}(x)A_{2}(x)+{{}^{\top}}B_{2}(x)B_{2}(x)\big{]}.$
■

We derive a new short proof of the next known theorem.

######  Theorem 4 ([14, Theorem E])

The following are equivalent:

  1. 1.

$S$ is a Stieltjes matrix moment sequence;

  2. 2.

$H_{m}$ and $EH_{m}$ are positive semi-definite for every $m\in\mathbb{N}$.

□

###### Proof

It suffices to show the converse sense. To this aim assume that $H_{m}$ and
$EH_{m}$ are positive semi-definite for every $m\in\mathbb{N}$, and let us
show that $S$ is a Stieltjes matrix moment sequence. Again from Theorem 1,
this is equivalent to show that $\mathbf{L}_{S}$ is a tracial moment
functional on $[0,+\infty[$, according to Proposition 1, this is equivalent to
show that $\mathbf{L}_{S}(P)\geq 0$ for all $P\geq 0$ on $[0,+\infty[$. Using
Lemma 2 it is enough to show that: $\mathbf{L}_{S}({{}^{\top}}PP)\geq 0$ and
$\mathbf{L}_{S}(X{{}^{\top}}PP)\geq 0$ for all
$P\in\mathcal{M}_{p}(\mathbb{R}[X])$. For
$P(X)=\displaystyle\sum_{i=0}^{n}A_{i}X^{i}\in\mathcal{M}_{p}(\mathbb{R}[X])$,
we have

${{}^{\top}}PP(X)=\displaystyle\sum_{i,j=0}^{n}{{}^{\top}}A_{i}X^{i+j}A_{j}.$

Since
$Tr({{}^{\top}}A_{i}S_{i+j+\epsilon}A_{j})=\displaystyle\sum_{k=0}^{m}{{}^{\top}}(A_{i}e_{k})S_{i+j+\epsilon}A_{j}e_{k},\quad
for\quad\epsilon\in\\{0,1\\}$, we get

$\mathbf{L}_{S}({{}^{\top}}PP)=\displaystyle\sum_{i,j=0}^{m}Tr({{}^{\top}}A_{i}S_{i+j}A_{j})=\displaystyle\sum_{k=0}^{m}\Big{(}\sum_{i,j=0}^{m}{{}^{\top}}(A_{i}e_{k})S_{i+j}A_{j}e_{k}\Big{)}\geq
0,$

and

$\mathbf{L}_{S}(X{{}^{\top}}PP)=\displaystyle\sum_{i,j=0}^{m}Tr({{}^{\top}}A_{i}S_{i+j+1}A_{j})=\displaystyle\sum_{k=0}^{m}\Big{(}\sum_{i,j=0}^{m}{{}^{\top}}(A_{i}e_{k})S_{i+j+1}A_{j}e_{k}\Big{)}\geq
0.$

■

###  3.3 Hausdorff matrix moment sequences.

The basic theorem in the description of non negative matrix polynomials on
$\left[0,1\right]$ states the following.

######  Theorem 5 ([7, Theorem 2.5])

The following are equivalent:

  1. 1.

$P\geq 0$ on $[0,1]$;

  2. 2.

    1. (a)

If $deg(P)$ is even, then
$P(X)={{}^{\top}}A(X)A(X)+X(1-X){{}^{\top}}D(X)D(X)$,

    2. (b)

If $deg(P)$ is odd, then $P(X)=X{{}^{\top}}A(X)A(X)+(1-X){{}^{\top}}B(X)B(X)$,

where $A$, $B$ are polynomials in $\mathcal{M}_{p}(\mathbb{R}[X]).$

□

As a consequence of the above characterisation of non negative polynomials on
$[0,1]$, the next results has been obtained in [7].

######  Theorem 6 ([7, Theorem 2.6])

Under the notation above. The following are equivalent:

  1. 1.

A family $(S_{0},\cdots,S_{n})$ is a Hausdorff moment sequence;

  2. 2.

    1. (a)

If $n=2m$, then $H_{m}$ and $(E-E^{2})H_{m}$ are positive semi-definite,

    2. (b)

If $n=2m+1$, then $EH_{m}$ and $(I-E)H_{m}$ are positive semi-definite.

□

Thanks to the expressions $X=X^{2}+X(1-X)$ and $1-X=(1-X)^{2}+X(1-X)$, we
derive the following finer characterisation of non negative matrix polynomials
improving Theorem 5.

###### Proposition 2

The following are equivalent:

  1. 1.

$P\geq 0$ on $[0,1]$;

  2. 2.

$P(X)=Q_{1}(X)+X(1-X)Q_{2}(X)$, where $Q_{1}$ and $Q_{2}$ are sums of
symmetric squares of polynomials in $\mathcal{M}_{p}(\mathbb{R}[X]).$

□

We deduce a short proof of the next improvement of Theorem 6.

###### Theorem 7

The following are equivalent:

  1. 1.

$S$ is a Hausdorff matrix moment sequence;

  2. 2.

$H_{m}$ and $(E-E^{2})H_{m}$ are positive semi-definite for every $m$.

□

###### Proof

As before, we assume that $H_{m}$ and $(E-E^{2})H_{m}$ are positive semi-
definite for every $m$, or equivalently that

$\displaystyle\sum_{0\leq i,j\leq m}c_{i}S_{i+j}{{}^{\top}}c_{j}\geq 0\mbox{
and}\displaystyle\sum_{0\leq i,j\leq
m}c_{i}(S_{i+j+1}-S_{i+j+2}){{}^{\top}}c_{j}\geq 0\mbox{ for
}c_{0},c_{1},...,c_{m}\in\mathbb{R}^{p},$

and let us show that $S$ is a Hausdorff matrix moment sequence. Again,
according to Proposition 1, this is equivalent to $\mathbf{L}_{S}$ is a
Tracial moment functional on $[0,1]$, from Theorem 1 this is equivalent to
show that $\mathbf{L}_{S}(P)\geq 0$ for all $P\geq 0$ on $[0,1]$. Now, using
Theorem 2, it is enough to show that $\mathbf{L}_{S}({{}^{\top}}PP)\geq 0$ and
$\mathbf{L}_{S}(X(1-X){{}^{\top}}PP)\geq 0$ for all
$P\in\mathcal{M}_{p}(\mathbb{R}[X])$. To this aim, for
$P(X)=\displaystyle\sum_{i=0}^{n}A_{i}X^{i}$, we have
${{}^{\top}}PP(X)=\displaystyle\sum_{i,j=0}^{n}{{}^{\top}}A_{i}X^{i+j}A_{j},$
and then

$\begin{array}[]{rl}\mathbf{L}_{S}({{}^{\top}}PP)&=\displaystyle\sum_{i,j=0}^{m}Tr({{}^{\top}}A_{i}S_{i+j}A_{j})=\displaystyle\sum_{k=0}^{m}\left(\sum_{i,j=0}^{m}{{}^{\top}}(A_{i}e_{k})S_{i+j}A_{j}e_{k}\right)\geq
0,\\\
\mathbf{L}_{S}(X(1-X){{}^{\top}}PP)&=\displaystyle\sum_{i,j=0}^{m}Tr\left({{}^{\top}}A_{i}(S_{i+j+1}-S_{i+j+2})A_{j}\right)\\\
\nobreak\leavevmode\hfil&=\displaystyle\sum_{k=0}^{m}\left(\sum_{i,j=0}^{m}{{}^{\top}}(A_{i}e_{k})(S_{i+j+1}-S_{i+j+2})A_{j}e_{k}\right)\geq
0.\end{array}$

■

##  4 Sequence of matrices defined by a linear recurrence relation

In this section, we are concerned with sequences of matrices defined by a
given $r$ initial data and satisfying a recursive linear relation of order
$r$. To be more precise, let $a_{0},a_{1},...,a_{r-1}$ be some fixed real
numbers with $a_{r-1}\neq 0$ and consider the sequence of matrices
$S=(S_{n})_{n\in\mathbb{N}}$ of order $p$ defined by the following linear
recurrence relation of order $r$,

$S_{n+1}=a_{0}S_{n}+a_{1}S_{n-1}+...+a_{r-1}S_{n-r+1},\quad n\geq r-1,$ (5)

where $S_{0},S_{1},...,S_{r-1}\in\mathcal{M}_{p}(\mathbb{R})$ are specified by
the initial conditions.

The case of scalar sequences, $p=1$, has been widely studied in the
literature, generally called $r$-generated Fibonacci sequences have been the
key point in solving several known problems; see for example [9] and [1].

Writing $S_{n}=(s_{i,j}^{(n)})_{1\leq i,j\leq p}$ for $n\in\mathbb{N}$, yields
(5) is equivalent to the following family of linear scalar relations,

$s_{i,j}^{(n+1)}=a_{0}s_{i,j}^{(n)}+a_{1}s_{i,j}^{(n-1)}+...+a_{r-1}s_{i,j}^{(n-r+1)},\quad
n\geq r-1,$ (6)

for every $1\leq i,j\leq p$.

The polynomial $Q_{S}(X)=X^{r}-a_{0}X^{r-1}-...-a_{r-2}X-a_{r-1}$ is called a
characteristic polynomial associated with the sequence
$S=(S_{n})_{n\in\mathbb{N}}$.

As in the scalar case, the characteristic polynomial is not unique.

###### Example 1

Let $M\in\mathcal{M}_{p}(\mathbb{R})$, and $S_{n}=(n+1)M$. Then the sequence
$S=(S_{n})_{n\in\mathbb{N}}$ verifies the two following recurrence relations:

$\begin{array}[]{lll}S_{0}=M,S_{1}=2M,S_{2}=3M&\mbox{ and
}\;S_{n+1}=S_{n}+S_{n-1}-S_{n-2},&\mbox{ for }n\geq 2,\\\
S_{0}=M,S_{1}=2M&\mbox{ and }\;S_{n+1}=2S_{n}-S_{n-1},&\mbox{ for }n\geq
1.\end{array}$

It follows that $P(X)=X^{3}-X^{2}-X+1$ and $Q(X)=X^{2}-2X+1$ are
characteristic polynomials for the sequence $S=(S_{n})_{n\in\mathbb{N}}$. □

Let $\mathcal{P}_{S}$ be the set of all characteristic polynomial associated
to the sequence $S=(S_{n})_{n\in\mathbb{N}}$.

###### Proposition 3

For every sequence $S=(S_{n})_{n\in\mathbb{N}}$ given by (5), their exists a
unique monic characteristic polynomial $P_{S}\in\mathcal{P}_{S}$ with minimal
degree. Moreover, every $P\in\mathcal{P}_{S}$ is a multiple of $P_{S}$. □

###### Proof

Let $W$ be the shift operator on
$\mathcal{L}\left(\mathcal{M}_{p}(\mathbb{R})^{\mathbb{N}}\right)$ defined by,

$W((S_{n})_{n\in\mathbb{N}})=(S_{n+1})_{n\in\mathbb{N}}.$

For $P(X)=\displaystyle\sum_{i=0}^{m}a_{i}X^{i}\in\mathbb{R}[X]$, we write
$P(W)=\displaystyle\sum_{i=0}^{m}a_{i}W^{i}$. Hence, it follows from (5) that

$Q\in\mathcal{P}_{S}\iff Q(W)(S)=0_{\mathcal{M}_{p}(\mathbb{R})}.$

We pose $I_{S}:=\left\\{P\in\mathbb{R}[X]:\quad
P(W)(S)=0_{\mathcal{M}_{p}(\mathbb{R})}\right\\},$ it is clear that $I_{S}$ is
an ideal in $\mathbb{R}[X]$ that is principal. Then there exists a unique
monic polynomial $P_{S}\in I_{S}$ of minimal degree such that

$\mathcal{P}_{S}=I_{S}=P_{S}\mathbb{R}[X].$

■

$P_{S}$ will be called the minimal polynomial associated with
$S=(S_{n})_{n\in\mathbb{N}}$.

###  4.1 Minimal Polynomial and Representing Matrix Measure

We characterize below the minimal polynomial associated with sequences given
by the recursive relation (5) in function to the minimal polynomials
associated with the entries. More precisely,

###### Theorem 8

Under the notation above, let
$S=(S_{n})_{n\in\mathbb{N}}=\big{(}(s_{i,j}^{(n)})_{1\leq i,j\leq
p}\big{)}_{n\in\mathbb{N}}$ satisfying (5) and $P_{ij}$ be the minimal
polynomial associated with $s_{ij}=(s_{i,j}^{(n)})_{n\in\mathbb{N}}$ for
$1\leq i,j\leq p$. Then

$P_{S}=lcm\big{\\{}P_{ij}|\quad 1\leq i,j\leq p\big{\\}}.$

where $lcm\big{\\{}P_{ij}|\quad 1\leq i,j\leq p\big{\\}}$ stands for the
lowest common multiple of the family of polynomials $\big{\\{}P_{ij}|\quad
1\leq i,j\leq p\big{\\}}$. □

###### Proof

Let $P_{ij}$ be the minimal polynomial associated with
$s_{ij}=(s_{i,j}^{(n)})_{n\in\mathbb{N}}$ for $1\leq i,j\leq p$, and
$P_{lmc}=lcm\big{\\{}P_{ij}|\quad 1\leq i,j\leq p\big{\\}}$. Let us show that
$P_{S}=P_{lmc}$. Since $P_{S}$ is a characteristic polynomial for
$s_{ij}=(s_{i,j}^{(n)})_{n\in\mathbb{N}}$, it is a multiple of $P_{ij}$ for
every $1\leq i,j\leq p$ and then $P_{S}$ is a multiple of $P_{lmc}$.
Conversely $P_{lmc}$ is a multiple of $P_{ij}$ for every $1\leq i,j\leq p$ and
then is a characteristic polynomial of
$s_{ij}=(s_{i,j}^{(n)})_{n\in\mathbb{N}}$ for every $1\leq i,j\leq p$. Now
since addition and scalar multiplication act linearly on entries, it follows
that $P_{lmc}$ is a characteristic polynomial for $S$. We deduce that
$P_{lmc}$ is a multiple of $P_{S}$. Which ends the proof. ■

We also have the next lemma. The proof is immediate and is left to the reader.

###### Lemma 3

Under the notation above, denote $Z(P)$ for the set of zeros of the polynomial
$P$. Then

$Z(P_{S}):=\displaystyle\bigcup_{1\leq i,j\leq p}Z(P_{ij}).$ (7)

□

We say that $P_{S}$ has distinct roots if
$Z(P_{S})=\big{\\{}\lambda_{1},\lambda_{2},...,\lambda_{k}\big{\\}}$, with
only simple roots.

###### Theorem 9

Let $S=(S_{n})_{n\in\mathbb{N}}$ be a recursive sequence, with
$Z(P_{S})=\left\\{\lambda_{1},\lambda_{2},...,\lambda_{k}\right\\}$. The
following are equivalent:

  1. 1.

$P_{S}$ has distinct roots;

  2. 2.

$S=(S_{n})_{n\in\mathbb{N}}$ admits a generating matrix-valued measure
$\sigma=(\sigma_{ij})_{1\leq i,j\leq p}$, with
$\sigma=\displaystyle\sum_{i=0}^{k}T_{i}\delta_{\lambda_{i}},$ and where
$T_{0},T_{1},...,T_{k}\in M_{p}(\mathbb{R}).$

□

###### Proof

If $S=(S_{n})_{n\in\mathbb{N}}$ admits a generating matrix-valued measure
$\sigma=(\sigma_{ij})_{1\leq i,j\leq p}$, then

$S_{n}=\big{(}(s_{i,j}^{(n)})_{1\leq i,j\leq
p}\big{)}_{n\in\mathbb{N}}=\displaystyle\int_{Supp(\sigma)}x^{n}d\sigma(x):=\left(\int_{Supp(\sigma)}x^{n}d\sigma_{ij}(x)\right)_{1\leq
i,j\leq p}..$

Which is equivalent to

$s_{ij}^{(n)}=\int_{Supp(\sigma)}x^{n}d\sigma_{ij}(x),\quad 1\leq i,j\leq
p,\quad n\in\mathbb{N}.$

Then the sequence $s_{ij}=(s_{ij}^{(n)})_{n\in\mathbb{N}}$ admits a generating
measure $\sigma_{ij}$, so according to (6) and the Proposition 2.4 in [1],
$P_{ij}$ has distinct roots. It follows that $P_{S}=P_{lcm}$ has also distinct
roots. Moreover, $Supp(\sigma_{ij})=Z(P_{ij}),$ then
$\sigma=(\sigma_{ij})_{1\leq i,j\leq p}$ is a representing measure for
$S=(S_{n})_{n\in\mathbb{N}}$ and by definition of $Supp(\sigma)$, we have

$Supp(\sigma):=\displaystyle\bigcup_{1\leq i,j\leq
p}Supp(\sigma_{ij})=\bigcup_{1\leq i,j\leq p}Z(P_{ij}).$ (8)

Thus from to (7) and (8), we obtain
$Supp(\sigma)=Z(P_{S})=\big{\\{}\lambda_{0},\lambda_{1},...,\lambda_{k}\big{\\}}$.
■

###### Example 2

  1. 1.

In Example 1, $S_{n}=(n+1)M$, the minimal polynomial $P_{S}(X)=(X-1)^{2}$ has
not distinct roots. Thus the sequence $S=(S_{n})_{n\in\mathbb{N}}$ admits no
representing measure.

  2. 2.

Let $A=\begin{pmatrix}5&-3\\\ 6&-4\end{pmatrix}$, we pose $S_{n}=A^{n}$ for
$n\in\mathbb{N}$, the characteristic polynomial of $A$ is given by;
$P_{A}(X)=(X+1)(X-2)=X^{2}-X-2$, so according to the Cayley-Hamilton theorem,
we have $S_{2}=S_{1}+2S_{0}$. Now setting $S_{0}=I_{2}$, $S_{1}=A$ and

$S_{n}=S_{n-1}+2S_{n-2},\quad n\geq 2,$ (9)

a direct computation using the Binet formula associated with recurrence
relation (9) on entries gives,

$S_{n}=\begin{pmatrix}-(-1)^{n}+2(2)^{n}&(-1)^{n}-(2)^{n}\\\
-2(-1)^{n}+2(2)^{n}&2(-1)^{n}-(2)^{n}\end{pmatrix}.$

Thus

$s_{ij}^{(n)}=\displaystyle\int t^{n}d\sigma_{ij}(t)\mbox{ for }1\leq i,j\leq
2.$

with
$\,\sigma_{11}=-\delta_{-1}+2\delta_{2},\,\sigma_{12}=\delta_{-1}-\delta_{2},\,\sigma_{21}=-2\delta_{-1}+2\delta_{2}$
and $\sigma_{22}=2\delta_{-1}-\delta_{2}.$

Finally $S_{n}=\displaystyle\int_{Supp(\sigma)}t^{n}d\sigma(t)$, with

$\sigma=\begin{pmatrix}\sigma_{11}&\sigma_{12}\\\
\sigma_{21}&\sigma_{22}\end{pmatrix}=\begin{pmatrix}-1&1\\\
-2&2\end{pmatrix}\delta_{-1}+\begin{pmatrix}2&-1\\\
2&-1\end{pmatrix}\delta_{2}$

and $Supp(\sigma)=Z(P_{A})=\\{-1,2\\}.$

□

In the remaining part, we consider recursive sequences whose minimal
polynomial have distinct roots and then admit a representing measure $\sigma$.
Our objective to find necessary and sufficient conditions so that $\sigma$ is
positive semi-definite representing finite atomic measure. In other words
$\sigma=\displaystyle\sum_{i=0}^{k}T_{i}\delta_{\lambda_{i}},$ with $T_{i}\geq
0$ for $i=0,...,k$.

Let $S=(S_{n})_{n\in\mathbb{N}}$ be a sequence of positive symmetric matrices
of order $p$ defined by the following linear recurrence relation (5) of order
$r$ . For $n\in\mathbb{N}$ we consider the Hankel $H(n)$ defined by:

$H(n)=\begin{pmatrix}S_{0}&S_{1}&\cdots&S_{n}\\\ S_{1}&S_{2}&\cdots&S_{n+1}\\\
\vdots&\vdots&\ddots&\vdots\\\ S_{n}&S_{n+1}&\cdots&S_{2n}\end{pmatrix}.$

###### Lemma 4

If $S=(S_{n})_{n\in\mathbb{N}}$ admits a positive semi-definite representing
finite atomic measure $\sigma$, then $H(r-1)\geq 0$. □

###### Proof

Suppose that $S=(S_{n})_{n\in\mathbb{N}}$ admits a positive semi-definite
representing finite atomic measure
$\sigma=\displaystyle\sum_{q=0}^{k}T_{q}\delta_{\lambda_{q}},$ with $T_{q}\geq
0$ for $q=0,...,k$.
Hence $S_{n}=\displaystyle\sum_{q=1}^{k}T_{q}\lambda_{q}^{n}$ for
$n\in\mathbb{N}$. Let $c_{0},...,c_{r-1}\in\mathbb{R}^{p}$, we have

$\displaystyle\displaystyle\sum_{i,j=0}^{r-1}c_{i}S_{i+j}{{}^{\top}}c_{j}$
$\displaystyle=$
$\displaystyle\sum_{i,j=0}^{r-1}c_{i}\left(\sum_{q=0}^{k}T_{q}\lambda_{q}^{i+j}\right){{}^{\top}}c_{j}$
$\displaystyle=$
$\displaystyle\sum_{q=0}^{k}\big{(}\sum_{i=0}^{r-1}c_{i}\lambda_{q}^{i}\big{)}T_{q}\big{(}\sum_{j=0}^{r-1}{{}^{\top}}c_{j}\lambda_{q}^{i}\big{)}$
$\displaystyle=$
$\displaystyle\sum_{q=0}^{k}P(\lambda_{q})T_{q}{{}^{\top}}P(\lambda_{q})\geq
0.$

With $P(X)=\displaystyle\sum_{i=0}^{r-1}c_{i}X^{i}$ denotes a vector
polynomial with coefficients $c_{0},...,c_{r}\in\mathbb{R}^{p}$. ■

We also have

###### Proposition 4

Let $S=(S_{k})_{k\in\mathbb{N}}$ be a sequence of positive matrices satisfying
(5) and let $\sigma=\displaystyle\sum_{q=0}^{k}T_{q}\delta_{\lambda_{q}}$ be a
generating matrix-valued measure of $S$. The following are equivalent:

  1. 1.

$H(r-1)\geq 0$;

  2. 2.

$H(n)\geq 0$ for all $n\geq r$.

□

###### Proof

We only show $(i)\Rightarrow(ii)$. Assume that $H(r-1)\geq 0$ and let
$c_{0},c_{1},...,c_{n}\in\mathbb{R}^{p}$ for $n\geq r$. We have

$\displaystyle\sum_{i,j=0}^{n}c_{i}S_{i+j}{{}^{\top}}c_{j}=\sum_{q=0}^{k}P(\lambda_{q})T_{q}{{}^{\top}}P(\lambda_{q}),$
(10)

where $P(X)=\displaystyle\sum_{k=0}^{n}c_{k}X^{k}$ denotes a vector polynomial
with coefficients $c_{0},...,c_{n}\in\mathbb{R}^{p}$. If we set
$c_{k}=(c_{i}^{(k)})_{1\leq i\leq p}$ for $0\leq k\leq n$, then
$P(X)=\big{(}P_{i}(X)\big{)}_{1\leq i\leq p}$ with
$P_{i}(X)=\displaystyle\sum_{k=0}^{n}c_{i}^{(k)}X^{k}\in\mathbb{R}_{n}[X].$
Now, for fixed $1\leq i\leq p$, there exists two polynomials $Q_{i}$ and
$R_{i}$ such that

$P_{i}=P_{S}Q_{i}+R_{i},\quad deg(R_{i})<deg(P_{S})=r.$

Hence

$P=P_{S}Q+R,$

with $Q(X)=\big{(}Q_{i}(X)\big{)}_{1\leq i\leq p}$ and
$R(X)=\big{(}R_{i}(X)\big{)}_{1\leq i\leq p}$. In particular
$P(\lambda_{q})=R(\lambda_{q})$ because $P_{S}(\lambda_{q})=0$ for $1\leq
q\leq k$. If we denote by
$\alpha_{0},\alpha_{1},...,\alpha_{r-1}\in\mathbb{R}^{p}$ the coefficients of
$R$ then according to the hypothesis and (10) we have

$\displaystyle\sum_{i,j=0}^{n}c_{i}S_{i+j}{{}^{\top}}c_{j}=\sum_{q=0}^{k}R(\lambda_{q})T_{q}{{}^{\top}}R(\lambda_{q})=\sum_{i,j=0}^{r-1}\alpha_{i}S_{i+j}{{}^{\top}}\alpha_{j}\geq
0.$

Finally $H(n)\geq 0.$ ■

###### Lemma 5

Let $S=(S_{n})_{n\in\mathbb{N}}$ be a sequence (5) of order $r$ and $P_{S}$ be
its minimal polynomial. If $H(r-1)\geq 0$, then

  1. 1.

The following are equivalent

    1. (a)

$P\in\mathcal{P}_{S}$;

    2. (b)

There exists $n\in\mathbb{N}^{*},P^{n}\in\mathcal{P}_{S};$

    3. (c)

For all integer $n\in\mathbb{N}^{*},P^{n}\in\mathcal{P}_{S}.$

  2. 2.

The minimal polynomial $P_{S}$ has distinct roots.

□

###### Proof

$1.$ $(a)\iff(c)$ Suppose that $H(r-1)\geq 0$. For
$P(X)=a_{m}X^{m}+\ldots+a_{1}X+a_{0}\in\mathbb{R}_{m}[X]$, $m\geq r$. We note
$\widehat{P}={{}^{\top}}(a_{0},\ldots,a_{m},0,\ldots)\otimes I_{p},$ where
$I_{p}$ is the identity matrix of order $p$.

$H(\infty):=(S_{i+j})_{i,j\geq
0}=\left(\begin{matrix}S_{0}&S_{1}&S_{2}&\ldots\\\ S_{1}&S_{2}&S_{3}&\ldots\\\
S_{2}&S_{3}&S_{4}&\ldots\\\ \vdots&\vdots&\vdots&\ddots\end{matrix}\right).$

From Lemma 4 we derive that $H(\infty)\geq 0$. Also we clearly have
${{}^{\top}}\widehat{X^{k}}H(\infty)\widehat{P}=a_{m}S_{k+m}+\ldots+a_{0}S_{k},$
$k\in\mathbb{N},$ because of (5) we obtain

$P\in\mathcal{P}_{S}\iff H(\infty)\widehat{P}=\widehat{0}.$ (11)

For $n\geq 2$, we show that $H(\infty)\widehat{P^{n}}=\widehat{0}\Rightarrow
H(\infty)\widehat{P^{n-1}}=\widehat{0}$. Indeed if
$H(\infty)\widehat{P^{n}}=\widehat{0}$, we get

$0_{\mathcal{M}_{p}(\mathbb{R})}={{}^{\top}}\widehat{P^{n-2}}H(\infty)\widehat{P^{n}}=H(\infty)\widehat{P^{2n-2}}=^{\top}\widehat{P^{n-1}}H(\infty)\widehat{P^{n-1}}.$

Thus $H(\infty)\widehat{P^{n}}=0\Rightarrow
H(\infty)\widehat{P^{n-1}}=\widehat{0}\Rightarrow...\Rightarrow
H(\infty)\widehat{P}=\widehat{0}.$ According to (11), $P\in\mathcal{P}_{S}$.
The converse direction follows from the fact that $\mathcal{P}_{S}$ is an
ideal.
$(c)\Rightarrow(b)$ Clear.
$(b)\Rightarrow(a)$ It is enough to use the fact that
$H(\infty)\widehat{P^{n}}=\widehat{0}\Rightarrow
H(\infty)\widehat{P^{n-1}}=\widehat{0}$.
$2.$ Let $P_{S}(X)=\displaystyle\prod_{i=0}^{r-1}(X-\lambda_{i})^{m_{i}}$ and
let $m=\displaystyle\max\\{m_{i}\mid i=1,\ldots,r\\}$. Then
$\displaystyle\left(\prod_{i=0}^{r-1}(X-\lambda_{i})\right)^{m}\in\mathcal{P}_{S}$,
because $\mathcal{P}_{S}$ is an ideal. Applying 1 of Lemma 5 , we obtain
$\displaystyle\prod_{i=0}^{r-1}(X-\lambda_{i})\in\mathcal{P}_{S}$. Since
$P_{S}$ is minimal, then
$P_{S}(X)=\displaystyle\prod_{i=0}^{r-1}(X-\lambda_{i})$, as desired. ■

We derive the following theorem.

###### Theorem 10

Let $S=(S_{n})_{n\in\mathbb{N}}$ be a sequence (5) of order $r$ and $P_{S}$ be
its minimal polynomial. The following assertions are equivalent:

  1. 1.

$P_{S}$ has distinct roots
$\left(P_{S}(X)=\displaystyle\prod_{k=0}^{r-1}(X-\lambda_{k})\right)$ and $S$
admits a positive semi-definite representing finite atomic measure
$\sigma=\displaystyle\sum_{i=0}^{r-1}T_{i}\delta_{\lambda_{i}};$

  2. 2.

$H(r-1)\geq 0.$

□

###### Proof

According to the Lemma 4 the direct sense is clear.
Now we assume that $H(r-1)\geq 0$, then by Lemma 5 the minimal polynomial
$P_{S}$ has distinct roots i.e
$P_{S}(X)=\displaystyle\prod_{k=0}^{r-1}(X-\lambda_{k})$, so from Proposition
4, we get $S=(S_{n})_{n\in\mathbb{N}}$ admit a generating finite atomic
measure $\sigma=\displaystyle\sum_{i=0}^{r-1}T_{i}\delta_{\lambda_{i}},$ with
$T_{0},T_{1},...,T_{r-1}\in M_{p}(\mathbb{R}).$ Then it remains to show that
$T_{i}\geq 0$ for $\displaystyle i=0...r-1$. We have $H(r-1)\geq 0$, then

$\displaystyle\sum_{i,j=0}^{r-1}c_{i}S_{i+j}{{}^{\top}}c_{j}=\sum_{k=0}^{r-1}P(\lambda_{k})T_{k}{{}^{\top}}P(\lambda_{k})\geq
0,$ (12)

with $P(X)=\displaystyle\sum_{k=0}^{r-1}c_{k}X^{k}$, $c_{k}\in\mathbb{R}^{p}$.
For $0\leq i\leq r-1$, denote $\displaystyle
L_{i}(X)=\prod_{\begin{subarray}{c}0\leq j\leq r-1\\\ j\neq
i\end{subarray}}\frac{(X-\lambda_{j})}{(\lambda_{i}-\lambda_{j})}$ and
$Q_{i}(X)=v.L_{i}(X),$ where $v$ is an arbitrary vector in $\mathbb{R}^{p}$.
It’s clear that $Q_{i}(\lambda_{j})=v\delta_{ij}$ with $\delta_{ij}$ is the
Kronecker delta and applying $\eqref{pos}$ for $P=Q_{i}$ we get

$vT_{i}{{}^{\top}}v\geq 0.$

The desired result is obtained. ■

###  4.2 Some special cases

We investigate next two special cases of recursive sequences. Without any loss
of generality, we can take $S_{0}=I_{p}.$

####  4.2.1 Case r=2.

Let $S=(S_{n})_{n\in\mathbb{N}}$ be a sequence of positive symmetric matrices
of order $p\geq 1$, such that $S_{n+1}=a_{0}S_{n}+a_{1}S_{n-1}$ for $n\geq
r-1=1$, $S_{0}=I_{p}$ and $S_{1}$ are given. If the characteristic polynomial
$P_{S}(X)=(X-\lambda_{0})(X-\lambda_{1})$ such that $\lambda_{0}<\lambda_{1}$,
according to Theorem 9, $S=(S_{n})_{n\in\mathbb{N}}$ admits a generating
matrix-valued measure
$\sigma=\displaystyle\sum_{i=0}^{1}T_{i}\delta_{\lambda_{i}},$ with
$T_{1},T_{2}\in\mathcal{M}_{p}(\mathbb{R})$. In particular

$\left\\{\begin{array}[]{lllcl}T_{0}&+T_{1}&=&I_{p},\\\
\lambda_{0}T_{0}&+\lambda_{1}T_{1}&=&S_{1},\\\ \end{array}\right.$

and then,

$T_{0}=\displaystyle\frac{1}{\lambda_{1}-\lambda_{0}}(\lambda_{1}I_{2}-S_{1})\mbox{
and
}T_{1}=\displaystyle\frac{1}{\lambda_{1}-\lambda_{0}}(S_{1}-\lambda_{0}I_{2}).$

We derive the next result

###### Theorem 11

Let $S=(S_{n})_{n\in\mathbb{N}}$ be recursive of order 2 with
$P_{S}(X)=(X-\lambda_{0})(X-\lambda_{1})$ its characteristic polynomial such
that $\lambda_{0}<\lambda_{1}$. The following are equivalent:

  1. 1.

$\sigma$ is positive semi-definite;

  2. 2.

$\lambda_{0}I_{2}\leq S_{1}\leq\lambda_{1}I_{2}$.

□

We derive the next corollary extending the scalar case.

###### Corollary 1

For any given positive symmetric matrix $A$ of order $p$, it exists a matrix
moment sequence $S=(S_{n})_{n\in\mathbb{N}}$ such that $S_{0}=I_{p}$ and
$S_{1}=A$. □

###### Proof

It suffices to choose $\lambda_{0}<\lambda_{1}$ such that
$\lambda_{0}I_{2}\leq A\leq\lambda_{1}I_{2}$ and $S=(S_{n})_{n\in\mathbb{N}}$
such that

$S_{n+1}=(\lambda_{0}+\lambda_{1})S_{n}-\lambda_{0}\lambda_{1}S_{n-1},\quad
n\geq 1,$

with $S_{0}=I_{p}$ and $S_{1}=A$. ■

In the scalar case $p=1$, we recover Proposition 4.1 of [1] in the following
corollary.

###### Corollary 2

A sequence $s=(s_{n})_{n\in\mathbb{N}}$ such that $s_{0}=1$,
$s_{1}\in\mathbb{R}$ are given and $s_{n+1}=a_{0}s_{n}+a_{1}s_{n-1}$, then the
following are equivalent:

  1. 1.

$s$ is a moment sequence;

  2. 2.

$P_{s}(s_{1})\leq 0$.

With $P_{s}$ is the associated characteristic polynomial. □

####  4.2.2 Case r=3.

Below we investigate the recurrence matrix sequences of order $3$. Let
$S=(S_{n})_{n\in\mathbb{N}}$ be a sequence of positive symmetric matrices of
order $p\geq 1$, such that $S_{n+1}=a_{0}S_{n}+a_{1}S_{n-1}+a_{2}S_{n-2}$ for
$n\geq r-1=2$, where $S_{0}=I_{p}$, $S_{1}$ and $S_{2}$ are given. In this
case using $(I,S_{1},S_{2})$ is a moment sequence if and only if
$S_{2}-S_{1}^{2}\geq 0$ and then there are initial data with no moment
sequence extension. We exhibit below constrictive method to get the
representing measure. We choose
$P_{S}(X)=(X-\lambda_{0})(X-\lambda_{1})(X-\lambda_{2})$ with
$\lambda_{0}<\lambda_{1}<\lambda_{2}$ as minimal polynomial of $S$ and let
$\sigma=\displaystyle\sum_{i=0}^{2}T_{i}\delta_{\lambda_{i}}$ be the
associated generating matrix-valued measure obtained from Theorem 9. The
coefficients $T_{1},T_{2},T_{3}\in\mathcal{M}_{p}(\mathbb{R})$, satisfy the
next Vandermonde system

$\left\\{\begin{array}[]{lllcl}T_{0}&+T_{1}&+T_{2}&=&I_{p},\\\
\lambda_{0}T_{0}&+\lambda_{1}T_{1}&+\lambda_{2}T_{2}&=&S_{1},\\\
\lambda_{0}^{2}T_{0}&+\lambda_{1}^{2}T_{1}&+\lambda_{2}^{2}T_{2}&=&S_{2}.\end{array}\right.$

Thus $T_{0},T_{1}$ and $T_{2}$ are given by the expressions,

$\left\\{\begin{array}[]{lll}T_{0}&=&{\delta}.(S_{2}-S_{1}(\lambda_{1}+\lambda_{2})+\lambda_{1}\lambda_{2}I_{p}),\\\
T_{1}&=&{\delta}.(-S_{2}+S_{1}(\lambda_{0}+\lambda_{2})-\lambda_{0}\lambda_{2}I_{p}),\\\
T_{2}&=&{\delta}.(S_{2}-S_{1}(\lambda_{0}+\lambda_{1})+\lambda_{0}\lambda_{1}I_{p}),\end{array}\right.$

where
$\delta=\left((\lambda_{2}-\lambda_{1})(\lambda_{2}-\lambda_{0})(\lambda_{1}-\lambda_{0})\right)^{-1}>0$.
Hence, we get the following result.

###### Proposition 5

The following are equivalent:

  1. 1.

$\sigma$ is positive semi-definite;

  2. 2.

$\left\\{\begin{array}[]{l}S_{2}-(\lambda_{1}+\lambda_{2})S_{1}+\lambda_{1}\lambda_{2}I_{p}\geq
0,\\\ -S_{2}+(\lambda_{0}+\lambda_{2})S_{1}-\lambda_{0}\lambda_{2}I_{p}\geq
0,\\\ S_{2}-(\lambda_{0}+\lambda_{1})S_{1}+\lambda_{0}\lambda_{1}I_{p}\geq
0.\end{array}\right.$

□

We derive a necessary condition to get $\sigma$ is positive semi-definite
given by

$\lambda_{0}I_{p}\leq S_{1}\leq\lambda_{2}I_{p}\mbox{ and
}\lambda_{0}^{2}I_{p}\leq S_{2}\leq\lambda_{2}^{2}I_{p}.$

###### Example 3

Let $S=(S_{n})_{n\in\mathbb{N}}$ be a sequence of positive matrices of order
$p=3$, such that $S_{n+1}=6S_{n}-10S_{n-1}+4S_{n-2}$ for $n\geq 2$,
$S_{0}=I_{3}=\begin{pmatrix}1&0&0\\\ 0&1&0\\\ 0&0&1\end{pmatrix}$,
$S_{1}=\begin{pmatrix}2&-1&0\\\ -1&2&-1\\\ 0&-1&2\end{pmatrix}$ and
$S_{2}=\begin{pmatrix}5&-4&1\\\ -4&6&-4\\\ 1&-4&5\end{pmatrix}$ are given. The
associated minimal polynomial is given by:

$P_{S}(X)=X^{3}-6X^{2}+10X-4=(X-\lambda_{0})(X-\lambda_{1})(X-\lambda_{2}).$

With $\lambda_{0}=2-\sqrt{2}$, $\lambda_{1}=2$ and $\lambda_{2}=2+\sqrt{2}$.
We have

$S_{2}-(\lambda_{1}+\lambda_{2})S_{1}+\lambda_{1}\lambda_{2}I_{3}=\begin{pmatrix}1&\sqrt{2}&1\\\
\sqrt{2}&2&\sqrt{2}\\\ 1&\sqrt{2}&1\end{pmatrix}\geq 0,$
$-S_{2}+(\lambda_{0}+\lambda_{2})S_{1}-\lambda_{0}\lambda_{2}I_{3}=\begin{pmatrix}1&0&-1\\\
0&0&0\\\ -1&0&1\end{pmatrix}\geq 0,$

$S_{2}-(\lambda_{0}+\lambda_{1})S_{1}+\lambda_{0}\lambda_{1}I_{3}=\begin{pmatrix}1&-\sqrt{2}&1\\\
-\sqrt{2}&2&-\sqrt{2}\\\ 1&-\sqrt{2}&1\end{pmatrix}\geq 0.$

Therefore, according to the proposition 5, $S=(S_{n})_{n\in\mathbb{N}}$ admits
a non negative matrix semi-definite representative measure $\sigma$. In this
case

$\sigma=T_{0}\delta_{\lambda_{0}}+T_{1}\delta_{\lambda_{1}}+T_{2}\delta_{\lambda_{2}},$

$T_{0}=\displaystyle\frac{1}{4}\begin{pmatrix}1&\sqrt{2}&1\\\
\sqrt{2}&2&\sqrt{2}\\\ 1&\sqrt{2}&1\end{pmatrix}$,
$T_{1}=\displaystyle\frac{1}{2}\begin{pmatrix}1&0&-1\\\ 0&0&0\\\
-1&0&1\end{pmatrix}$ and
$T_{2}=\displaystyle\frac{1}{4}\begin{pmatrix}1&-\sqrt{2}&1\\\
-\sqrt{2}&2&-\sqrt{2}\\\ 1&-\sqrt{2}&1\end{pmatrix}$.

In particular
$S_{n}=T_{0}\lambda_{0}^{n}+T_{1}\lambda_{1}^{n}+T_{2}\lambda_{2}^{n}$ for all
$n\in\mathbb{N}.$ □

###### Remark 1

If $S=(S_{n})_{n\in\mathbb{N}}=\big{(}(s_{i,j}^{(n)})_{1\leq i,j\leq
p}\big{)}_{n\in\mathbb{N}}$ is a moment sequence associated with a non
negative measure $\sigma=(\sigma_{ij})_{1\leq i,j\leq p}$, then $s_{ii}^{(n)}$
is also a moment sequence. Indeed, we have $\sigma_{ii}\geq 0$ and
$s_{ii}^{(n)}=\displaystyle\int x^{n}d\sigma_{ii}(x)$ for $n\geq 0$.

Moreover Example 3 shows that this no more true for $s_{ij}$ when $i\neq j$.
On the other hand the next sequence

$S_{n}=\begin{pmatrix}1+2^{n}&1+3.2^{n}\\\ 1+3.2^{n}&2^{n}\\\
\end{pmatrix}=\begin{pmatrix}1&1\\\ 1&0\\\
\end{pmatrix}1^{n}+\begin{pmatrix}1&3\\\ 3&1\\\ \end{pmatrix}2^{n}$

has all its entries moment sequences but is not a moment sequence since
$\begin{pmatrix}1&3\\\ 3&1\\\ \end{pmatrix}$ and $\begin{pmatrix}1&1\\\ 1&0\\\
\end{pmatrix}$ are not positive. □

## References

  * [1] R. BenTaher, M. Rachidi, and E. H. Zerouali, Recursive subnormal completion and the truncated moment problem, Bulletin of the London Mathematical Society 33 (2001), no. 4, 425–432.
  * [2] C. H. Chidume, M. Rachidi, and E. H. Zerouali, Solving the general truncated moment problem by the $r$-generalized Fibonacci sequences method, Journal of mathematical analysis and applications 256 (2001), no. 2, 625–635.
  * [3] J. Cimprič and A. Zalar, Moment problems for operator polynomials, Journal of Mathematical Analysis and Applications 401 (2013), no. 1, 307–316.
  * [4] R. Curto and L. Fialkow, Solution of the truncated complex moment problem for flat data, in Memoirs of the American Mathematical Society, vol. 568, American Mathematical Society, Providence (1996).
  * [5] R. Curto and L. Fialkow, Flat extensions of positive moment matrices: recursively generated relations, in Memoirs of the American Mathematical Society, vol. 648, American Mathematical Society, Providence (1998).
  * [6] R. Curto and L. Fialkow, The truncated complex $\mathbf{K}$-moment problem, Transactions of the American mathematical society 352 (2000), no. 6, 2825–2855.
  * [7] H. Dette and W. J. Studden, Matrix measures, moment spaces and Favard’s theorem for the interval $[0,1]$ and $[0,+\infty[$, Linear Algebra and its Applications 345 (2000), no. 1-3, 169–193.
  * [8] M. A. Dritschel and J. Rovnyak, The operator Fejér-Riesz theorem, in A glimpse at Hilbert space operators, Theory Adv. Appl. 207 (2010), 223–254.
  * [9] F. Dubeau, W. Motta, M. Rachidi, and O. Saeki, On weighted $r$-generalized Fibonacci sequences, Fibonacci Quarterly 35 (1997), no. 2, 102–110.
  * [10] H. Hamburger, Über eine erweiterung des stieltjesschen momentenproblems, Mathematische Annalen 81 (1920), no. 2, 235–319.
  * [11] F. Hausdorff, Summationsmethoden und momentfolgen. I, Mathematische Zeitschrift 9 (1921), no. 1, 74–109.
  * [12] D. Kimsey and H. Woerdeman, The truncated matrix-valued $\mathbf{K}$-moment problem on $\mathbb{R}^{d}$, $\mathbb{C}^{d}$ and $\mathbb{T}^{d}$, Transactions of the American Mathematical Society 365 (2013), no. 10, 5393–5430.
  * [13] D. Kimsey and M. Trachana, On a Solution of the Multidimensional Truncated Matrix-Valued Moment Problem, Milan Journal of Mathematics 90 (2022), no. 1, 1–85.
  * [14] M. Krasnoselskiı and M. G. Krein, Fundamental theorems on the extension of Hermitian operators and certain of their applications to the theory of orthogonal polynomials and the problem of moments, Uspehi Mat. Nauk 2 (1947), no. 3(19), 60–106.
  * [15] M. G. Krein, Fundamental aspects of the representation theory of hermitian operators with deficiency index $(m,m)$, AMS Translations, series 97 (1970), no. 2, 75–143.
  * [16] C. T. Le, Tracial moment problems on hypercubes, Operators and Matrices 14 (2020), no. 4, 1015–-1027.
  * [17] M. Riesz, Sur le problème des moments, Arkiv för Mat., Astr. och Fysik 17 (1923), no. 12, 1–52.
  * [18] K. Schmüdgen The moment problem, The moment problem, Vol. 9, Berlin: Springer, 2017.
  * [19] T. J. Stieltjes, Recherches sur les fractions continues, Annales de la Faculté des sciences de Toulouse: Mathématiques 8 (1894), no. 4, J1–J122.
  * [20] J. Stochel, Solving the truncated moment problem solves the full moment problem, Glasgow Mathematical Journal 43 (2001), no. 3, 335–341.
  * [21] F.-H. Vasilescu, Dimensional stability in truncated moment problems, Journal of Mathematical Analysis and Applications 388 (2012), no. 1, 219–230.