Quantum SVR for Chlorophyll Concentration Estimation in Water With Remote Sensing

The increasing availability of quantum computers motivates researching their potential capabilities in enhancing the performance of data analysis algorithms. Similarly, as in other research communities, also in remote sensing (RS), it is not yet defined how its applications can benefit from the usage of quantum computing (QC). This letter proposes a formulation of the support vector regression (SVR) algorithm that can be executed by D-Wave quantum computers. Specifically, the SVR is mapped to a quadratic unconstrained binary optimization (QUBO) problem that is solved with quantum annealing (QA). The algorithm is tested on two different types of computing environments offered by D-Wave: the advantage system, which directly embeds the problem into the quantum processing unit (QPU), and a hybrid solver that employs both classical and QC resources. For the evaluation, we considered a biophysical variable estimation problem with RS data. The experimental results show that the proposed quantum SVR implementation can achieve comparable or, in some cases, better results than the classical implementation. This work is one of the first attempts to provide insight into how QA could be exploited and integrated in future RS workflows based on machine learning (ML) algorithms.

marketing, sociology, epidemiology, and risk analysis [4]. In the field of remote sensing (RS), the regression analysis has been applied in different applications [5], [6].
In the context of quantum machine learning (QML) [7], [8], [9], only few works have already addressed regression analysis problems. For instance, a quantum version of a linear regression algorithm and of a ridge regression algorithm have been proposed [10], [11]. Among the different paradigms of quantum computing (QC), quantum annealing (QA) has recently provided promising results in diverse ML applications [12], [13]. QA is a metaheuristic for solving combinatorial optimization problems [14], [15]. QA is closely related to adiabatic QC (AQC) [16], which was shown to be polynomially equivalent to the universal gate-based model, which is a different paradigm of QC [17]. However, QA can only solve a specific class of problems, and therefore, the redefinition of ML algorithms in a suitable format is one of the central design challenges when working with QA-enhanced ML models [18]. This work presents an implementation of the support vector regression (SVR) [19] algorithm that uses QA for solving the optimization problem related to the training phase of the SVR algorithm. Previous works tried to apply QA to optimize the training procedure of a support vector machine for classification tasks [13].
A similar implementation of QA-optimized SVR algorithm was proposed for facial landmarks detection [20]. Specifically, our implementation uses a similar workflow for constructing the quadratic unconstrained binary optimization (QUBO), but the mathematical formulation presents some differences in the constraints enforcement procedure (Section II). In addition, in this work, we propose six different methods to combine the solutions returned by the annealer when running the problem on the advantage system (Section II). Moreover, our implementation was tested on both hybrid and direct quantum processing unit (QPU) solvers, whereas [20] tested the quantum SVR (QSVR) only on hybrid solvers. This work focuses on biophysical parameter estimation related to chlorophyll concentration in water [21], [22], [23]. The proposed implementation was tested on a synthetic and a real RS dataset related to chlorophyll concentration in water. The quantum system used in the experiments was provided by the company D-Wave. Specifically, the experimental validation was conducted on the D-Wave Advantage_system4.1 solver and the hybrid_binary_quadratic_model_version2 hybrid solver. The access to such computational resources was provided through the D-Wave Leap cloud service. The purpose of this work 1558-0571 © 2022 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission.
is to investigate how QA could improve the existing ML frameworks for RS applications.

A. Support Vector Regression
The mathematical formulation of the -insentive SVR is now briefly described. Let T = {(x n , y n ), n = 0, . . . , N − 1} be the dataset used for the training phase constituted by N training samples. Each of such samples is formed by a feature vector x n ∈ R d , where d is the dimension of the feature space, and a target value y n ∈ R. It can be shown, with some mathematical manipulation, that the training phase amounts to the solving of the following constrained optimization problem: where the constraints to be enforced are In the optimization problem, the terms α α α = {α n : n = 0, . . . , N − 1} andα α α = {α n : n = 0, . . . , N − 1} are the variables with respect to which the problem is optimized. The terms C and are parameters related to the SVR that controls the overfitting and the error sensitivity, respectively. After finding the values of α α α andα α α, the prediction function is defined as The value of b can be deducted from any data point for which 0 < α n < C or 0 <α n < C using the formula To obtain a more robust estimation of b, it is preferable, however, to average the result from multiple data points [24].
The term k(x n , x m ) indicates the kernel function; in our experiments, a radial basis function (RBF) kernel, whose formula is given by: e (−γ ||x n −x m || 2 ) , has been employed.

B. QUBO Problem Formulation
To be processed by the quantum annealer, a problem must be in the form of either a Ising spin problem [25] or a QUBO. For our purposes, the original optimization problem related to the SVR was turned into a QUBO problem. Such a problem can be expressed according to the following formula: Terms a i ∈ {0, 1} are the binary variables of the QUBO problem, and Q is an upper-diagonal matrix called QUBO weight matrix that defines the problem. To turn the original problem into a QUBO, it is first necessary to encode the variables α α α andα α α in the binary variables a i . To do so, the following encoding strategy is applied [13]: In (6) and (7), the value K corresponds to the number of logical qubits used to encode each variable, whereas B is the value of the base used for the encoding. From the abovementioned equations, it is possible to note that the total number of variables of the QUBO problem is 2K N, and that the first K N variables are used to encode the α α α variables, whereas the last K N ones are used to encode the variablesα α α. The parameter P ≥ 0 is used for enabling the usage of negative exponents in the encoding procedure. To enforce the constraint defined in (2a), a square penalty term whose strength is regulated by the hyperparameter ξ is added to the cost function. The constraints in (2b) and (2c) are implicitly satisfied by the choice of the hyperparameter C: from the equations, it is possible to see that the maximum value that each variable can take is Therefore, by choosing a value of C equal or higher than this quantity, it is possible to guarantee the enforcing of the constraints. The lower bound is always satisfied, because each α n andα n are non-negative by definition. Moreover, another penalty term, whose influence is controlled by the hyperparameter β, is added to enforce that, for each value of n, at least one of α n orα n is equal to 0, or equivalently: By adding the penalty terms to the cost functions and by applying the encoding equations, it is possible to obtain the final formulation of the QUBO problem The termQ is a 2K N × 2K N matrix that defines the problem, whose elements are given bỹ with n, m ∈ {0, . . . , N − 1}, i, j ∈ {0, . . . , K − 1}, and s, t ∈ {0, 1}. Since the QUBO weight matrix Q is upper triangular, it is obtained fromQ using the formula This problem formulation is used for the QPU of the advantage system and the leap's hybrid solver, as they both input a binary quadratic problem. In the case of leap's hybrid solver, the problem is optimized using both classical and QC resources. The allocation of such resources and the problem decomposition is done automatically by the solver.
The final step for optimizing a problem with the QPU is the minor embedding [26]. In this process, the QUBO problem is embedded in the hardware architecture used for annealing. Specifically, each logical variable a i is mapped to a chain of qubits, i.e., a group of connected qubits used to represent a specific logical variable of the QUBO problem. The reason why chains are needed is that it is not always possible to directly map the optimization problem directly into the hardware topology. The values of the elements Q i, j of the Q matrix, which represent the coefficients associated with the terms a i a j of the cost function, are mapped to physical connections of chains of qubits.

C. Solutions Combination Techniques
The advantage system outputs 10 000 reads (i.e., solutions) with different levels of energy. In this work, we select the 40 solutions with the lowest energy, and we fuse them to compute the final solution (the number 40 is chosen arbitrarily and other options can be considered). To combine them, we propose six approaches based on different weighted average formulas. The predictions of the 40 solutions are then evaluated on the training dataset. To each of them is assigned a score depending on the value of a specific loss function, and such scores are then used to obtain the coefficients of the weighted average. Specifically, the scores are calculated by considering the multiplicative inverse of the value of a loss function between the actual and the predicted value; therefore, a lower value of the loss will be associated with a higher score value. The loss functions considered are even and nonnegative. The combinations methods differ for the choice of the loss function and how the scores are used to get the coefficients. Every method ensures that each weight coefficient is non-negative; its value is lower or equal than 1, and that their sum is equal to 1. A brief description of the used methods is now provided.

5) QSVR 5:
Only the best solution in terms of MSE is considered, and this is done by setting the weights associated with the best solution to 1 and all the others to 0. 6) QSVR 6: To each solution is assigned the same weight; therefore, a simple average is performed.

III. DATASET DESCRIPTION
In this section, a brief overview and description of the dataset used for the experimental analysis are provided.
1) MERIS: The first dataset used is a synthetic dataset whose aim is to simulate the concentration of chlorophyll concentration and its relation to optical measurements. The wavelengths considered are the first eight spectral bands of the Multispectral Medium Resolution Imaging Spectrometer (MERIS) sensor (412.5, 442.5, 490, 510, 560, 620, 665, and 681.25 nm). The procedure employed to generate the dataset is the one described in [27]. 2) SeaBAM: The dataset contains information about 919 measurements regarding chlorophyll-a water concentration performed in Europe and United States. The value of concentration ranges between 0.019 and 32.787 mg/m 3 . The sensor used for the measurements is the sea-viewing wide field of view sensor (SeaWIFS), and the wavelengths considered in the experiments were 412, 443, 490, 510, and 555 nm. In both cases, the feature vector is constructed by considering the spectral measures at different wavelengths, whereas the target value is the corresponding chlorophyll concentration.

IV. EXPERIMENTAL ANALYSIS
For each dataset, two experiments were conducted: one using the D-Wave advantage QPU and the other using the leap's hybrid binary quadratic model (BQM) solver. The implementation of the classical SVR was done using the Python library sci-kit learn. In each setting, ten test  II   VALUES OF MSE OBTAINED ON THE TEST SET FOR DIFFERENT COMBINATION METHODS AND RUNS ON THE SYNTHETIC DATASET OBTAINED BY THE  ADVANTAGE SYSTEM. VALUES IN BOLD INDICATE THE RESULTS WHERE THE QUANTUM IMPLEMENTATION PERFORMED BETTER  runs were carried out; each one using different datasets for training and testing that were randomly chosen from the initial dataset. In each problem instance, the results were compared with a traditional SVR on the same datasets. Moreover, the hyperparameters for the quantum and the classical implementation of the SVR were the same for each test run. For the experiments using the advantage solver, the number of training samples was 30, whereas for the hybrid solver, the number was 50. The reason for this is that the advantage system could not always find an embedding with bigger problem instances. This is likely due to the structure of the problem itself that presents many interactions between variables that makes difficult finding an embedding as the number of variables increases [28]. In each test run, the hyperparameter γ was validated classically using a validation dataset and an SVR, and the validation dataset was divided into two parts: the first one was used to train a classical SVR with a given hyperparameter configuration values, while the other was used for testing. The configuration that achieved the best performances in terms of MSE was then used in the test run. The values of γ were selected from the range [0.1, 0.5, 1, 1.5, 2, 3, 4, 5, 7, 10, 20, 50], whereas the value of C was set to C min , with C min being the quantity defined in (8), which is equal to the maximum value each α n and α n can take. The values related to the problem encoding were B = 4 and P = 1 for the synthetic dataset and B = 5 and P = 0 for the Sea-viewing Bio-optical Algorithm Mini-Workshop (SeaBAM), whereas the number of logical qubits K used to encode the original problem variables was equal to 2 in each test run. For the training phase with the SeaBAM dataset, the values of both the feature vector and the target value were converted to the logarithmic domain as was done in [21], [22], and [23]. The reason for this is that the distribution of the biophysical quantities is assumed to be log-normally distributed [29]. Table I reports the result in terms of MSE obtained by the hybrid solver on both the synthetic and the seaBaM dataset, while the results obtained by the advantage system for the synthetic and SeaBAM datasets are reported in Tables II and III, respectively. When considering the experiments on the synthetic datasets, the quantum SVR achieved similar results to its classical counterpart on both the hybrid solver and the advantage system. In the experiments on the SeaBAM dataset, the QSVR on average performed better than the classical one on the hybrid solver, whereas in the experiments using the advantage system, the classical implementation performed slightly better, but the quantum version managed to obtain good results nevertheless and to perform better on some test runs. The different versions of QSVR performed similarly in the synthetic dataset; while on the SeaBAM dataset, there was more variation within the results, and this might indicate that the correct choice for the solutions combination technique becomes more important, as the data complexity increases. In the experiments with the SeaBAM dataset, the QSVR 5 managed to obtain the best results in average among the different QSVR implementations, but in some specific problem instances, it performed worse compared with the other solutions combination techniques. The QSVR 2 obtained the second best average results among the quantum implementations, but it was the implementation that outperformed the classical one the most number of times, five out of ten. The link to the repository associated with this work can be found at. 1

V. CONCLUSION
The main objective of this work was to investigate how QA could enhance an SVR algorithm for an RS application. The proposed algorithm was tested on both the D-Wave Advantage system and on the hybrid solver. The results show that the quantum implementation of SVR could achieve similar or, in some cases, even better results than the classical SVR. This is indicative of the potential of QA especially when considering that the original problem was continuous and unconstrained, and that it had to be modified and adapted to be solved by the annealer. In general, the hybrid solver provided better results than the Advantage system, and it could also solve bigger problem instances. Therefore, in the near future, practical applications will likely run on a hybrid framework.