Journal of Computer Chemistry, Japan -International Edition
Online ISSN : 2189-048X
ISSN-L : 2189-048X
Lattice Folding Simulation of Peptide by Quantum Computation
Rui SAITO, Koji OKUWAKI, Yuji MOCHIZUKI, Ryutaro NAGAI, Takumi KATO, Kenji SUGISAKI, Yuichiro MINATO
Author information
JOURNAL OPEN ACCESS FULL-TEXT HTML
Supplementary material

2023 Volume 9 Article ID: 2022-0036

Details
Download PDF (2858K)
Download citation RIS

(compatible with EndNote, Reference Manager, ProCite, RefWorks)

BIB TEX

(compatible with BibDesk, LaTeX)

Text
How to download citation
Contact us
Abstract

Computational protein folding has attracted considerable interest over the years, including molecular simulations and artificial intelligence assisted methods. On the other hand, research and development of quantum computer hardware and software have been thriving recently. In this paper, we report a case study of peptide (PSVKMA) folding based on a two-dimensional lattice model, by using both the blueqat quantum simulator (called AutoQML) and the IonQ quantum device. As a result, it was found that the actual device was still susceptible to noises.

1 INTRODUCTION

About 40 years have passed since the proposal of the quantum computer by Feynman [1], and various actual machines as well as simulators have been developed [2], raising expectations for its great potential. In particular, quantum chemical computation, which was pioneered by Ref [3] is being actively researched and developed as a killer application for molecular designs and material developments [4, 5]. On the other hand, Aspuru-Guzik et al. have attempted to apply quantum computations to the folding problem of proteins within lattice models [6,7,8], from the interest of theoretical drug discovery and biophysics. There have been recent reports in this field as well [9,10,11,12]. Here, the interactions between amino acid residues are represented by a coarse-grained (CG) picture based on the Miyazawa-Jernigan (MJ) model [13, 14] or the hydrophobic-polar (HP) model [15], and their optimal arrangement on the lattice may be solved by a quantum annealing approach on annealer type devices just such as D-Wave [16].

The quadratic unconstrained binary optimization (QUBO) [17] is a well-known method for defining optimization problems in quantum annealing. QUBO was actually used in folding simulations [6,7,8, 12, 18], however it can have a drawback that it becomes cumbersome when representing many-body interactions [19]. This difficulty of the setup for mapping a given problem to an actual quantum circuit can increase rapidly. We thus employed a cloud environment for quantum computations named blueqat AutoQML [20]. This web-based environment provides conveniently automated conversions from QUBO to quantum approximate optimization algorithm (QAOA) [21, 22] with variational quantum eigensolver (VQE) [23] for gate type devices, by which the difficulty of lengthy and complicated setups for QUBO solving can be avoided. Furthermore, the blueqat AutoQML is directly connectable to actual quantum computer devices such as IonQ of ion-trap type [24].

In the precedent short article [25], we briefed a couple of trial folding simulations (of three and five qubits) for the PSVKMA peptide (Pro-Ser-Val-Lys-Met-Ala) [7, 8], by using both blueqat AutoQML and IonQ. This article reports details of the above two trials and also another trial of five qubits. A convenience of conversion from QUBO to QAOA will be demonstrated by the well-defined folding problem [7, 8].

2 SET-UP OF TRIAL SIMULATIONS

As mentioned above, the blueqat AutoQML [20] is a web-based application by which a given QUBO equation is automatically converted to an optimization problem of QAOA and the corresponding solution samples are obtained. A simple illustration of “q0q1” interaction term may be embedded to a quantum circuit of RZ and CNOT as shown in Figure 1. Although QUBO formally contain up to two-body terms in quantum annealing, the cost Hamiltonian could have many-body (or more than two-body) terms, and these higher terms should be reduced to two-body terms by introducing ancilla qubits, leading to more demands of qubits. In contrast, many-body terms can be handled without ancilla qubit in QAOA. The blueqat AutoQML supports many-body terms appeared in the cost Hamiltonian; n-body terms may be expressed with RZ and CNOT gates [26]. Once the solution of QAOA is obtained, the corresponding quantum circuit can be directly passed to actual quantum computers just as the present IonQ [24] in a cloud environment; we actually used an 11-qubit device composed of 13 171Yb+ ions [27]. As a whole, a user-friendly interface is provided by the blueqat AutoQML.

Figure 1.

 Quantum circuit for interaction of two qubits.

The qubit treatment of protein folding within a two-dimensional Cartesian lattice model [7, 8] is resumed as below. As illustrated in Figure 2, a certain folding is presented by a series of turns, and each turn is given by a qubit pair, where “00,” “01,” “10” and “11” correspond to “turn down,” “turn right,” “turn left” and “turn up,” respectively. A qubit sequence with fixing the first pair “01” may thus be written as    q = 010q1q2q3…q(2N-6)q(2N-5),    (1)    where “0q1” specifies the second turn and N is a formal number of residues (depending on the problem setting). The energy function to be minimized contains constraint terms as the following equation.    E(q) = EBack(q)+EOverlap(q)+EPair(q)    (2)    EBack(q) is the term that penalizes the amino acid residues if they backward overlap, and EOverlap(q) is a similar term for the forward overlap. These two terms lead to destabilizations to prevent unphysical folding. In contrast, EPair(q) presents stabilizations due to residue-residue interactions, where the pair parameter values for PSVKMA [7, 8] are P:K = −1, P:A = −2, S:M = −3, and V:A = −4 by MJ [14]. These inter-residue parameters could effectively describe several stabilizations of such as hydrogen bond and hydrophobic attraction, which are essential for protein folding.

Figure 2.

 Relation between turn and qubit pair.

Three trials of quantum folding simulation were conducted, as shown in Figure 3; Trial 1 and Trial 2 were briefed in Ref [25]. Note that there is no direct correspondence to the experimental scheme names and qubit sequences used in Ref [7]. due to several reasons. Figure 4 illustrates the final correct structures presumed, where the stabilization interactions by EPair(q) is indicated with multiple bar.

Figure 3.

 Three trial models for folding on 2-dimensional lattice (see text). Trial 1 has three qubits to be determined, and Trials 2 and 3 are of five qubits problem.

Figure 4.

 Presumed correct folding structures of three trial models. Multiple bar indicates stabilization of residue-residue interaction between two residues. The corresponding qubit sequences to be solved for Trial 1, Trial 2 and Trial 3 are “010,” “01010” and “01011,” respectively.

The formal N of equation (1) is 6 for the PSVKMA peptide, and the number of required qubits is then 7, leading to too lengthy and complicated expression for the free folding [7, 8]. By fixing a certain part of peptide, the problem complication and size can be reduced. For the simplest Trial 1 in which PSVK is already folded, three qubits (q0, q1 and q2) are to be determined for the remaining MA, and the corresponding energy function to be optimized is given with the condition of q∈{0,1} as follows.    E(q) = − 1 − 4q2 + 9q0q2 + 9q1q2 − 16q0q1q2    (3)    There is a product term of three qubits (or three-body interaction). The energy expression of Trial 2 with the pre-folded PSV part is given as equation (4).    E(q) = − 3q1 + 7q0q1 + 18q1q2 − 15q0q1q2 − 4q0q3 − 2q1q3 + 15q0q1q3 + 7q2q3 + 4q0q2q3 − 7q1q2q3 − 24q0q1q2q3 + 7q1q4 + 4q0q1q4 − 18q1q2q4 + 7q3q4 + 4q0q3q4 − 7q1q3q4 − 24q0q1q3q4 − 18q2q3q4 + 20q1q2q3q4 + 29q0q1q2q3q4    (4)    This five qubits expression is considerably long relative to equation (3) of three qubits, where up to five-body terms appear. Similarly, the expression of Trial 3 has the following form.    E(q) = 4q0q1 + 15q1q2 − 15q0q1q2 − 4q0q3 + 2q0q1q3 + 7q2q3 + 4q0q2q3 − 20q1q2q3 + 9q0q1q2q3 + 7q1q4 + 4q0q1q4 − 18q1q2q4 + 7q3q4 + 4q0q3q4 − 20q1q3q4 + 9q0q1q3q4 − 18q2q3q4 + 53q1q2q3q4 − 33q0q1q2q3q4    (5)   

The blueqat AutoQML [20] could accept binary optimization formulas of equations (3), (4) and (5) and convert the corresponding QAOA equations to obtain the list of probabilities of qubit solutions.

3 RESULTS AND DISCUSSION

Although the results of Trial 1 were already presented in the previous short article [25], they are here restated for the convenience of readers. Figure 5 captures a web interface image of the blueqat AutoQML for Trial 1, where equation (3) is visible in the “QUBO” window. The step size of 5 means the number of iteration steps in a QAOA circuit of the blueqat simulator, where larger number may improve the accuracy of solution at the cost of deeper circuit. The correct answer of Trial 1 is “q0q1q2 = 001” (see again Figure 4), and it was successfully evaluated with a high probability of 0.99. The quantum circuit was then available as partly shown in the bottom window of Figure 6. This information was put forward to the IonQ quantum device, and the results of run were obtained as a number of counts for the respective qubit sequences; a python script for IonQ is listed in Supporting Information (SI). Figure 7 captures the actual results of 100 shots in a single batch at IonQ. Table 1 compares the probabilities of qubit sequences. In the IonQ results, the probability of “010” is the highest (or 0.41), however other non-preferable sequences with destabilization penalty such as “011” and “101” have sizable probability values. suggesting the affect from noises (e.g., single qubit-flip).

Figure 5.

 Captured web interface image of blueqat AutoQML simulator for Trial 1 [25]. Both input QUBO equation and optimized results are visible (refer to Table 1 as well).

Figure 6.

 Captured web interface image of quantum circuit output from blueqat AutoQML simulator for Trial 1. This information can be passed to IonQ quantum device straightforwardly.

Figure 7.

 Captured web interface image of python script and result by first IonQ batch (100 shots) [25]. The run result is visible in “Counter” line (as ‘001’: 41, ‘111’: 13, ‘011’: 12, ‘000’: 10, ‘100’: 10, ‘110’: 6, ‘101’; 6, ‘010’: 2). Here, two long blank (gray color) are manually put to hide security keys for actual usage of IonQ; “Internet address” line and “api = bqcloud.register_api” line.

Table 1.  Comparison of probabilities for Trial 1 [25].a)
Seq. blueqat IonQ
001 0.991 0.41
111 0.004 0.13
000 0.003 0.10
110 0.002 0.06
101 0.0 0.06
011 0.0 0.12
100 0.0 0.10
010 0.0 0.02

a)The number of shots in the first IonQ run was 100, and thus the normalization is simple.

For Trial 2, the blueqat AutoQML simulations with the step size of 5 were done five times, and the results are given in SI because of bulky amounts. Although the correct qubit sequence of “01010” was observed as the highest probability (its maximum value was 0.345), other unfavorable qubit sequences had probability values from 0.1 to 0.2. The simulations with the step size of 1 were performed ten times, and the corresponding results varied considerably, as listed in SI. Nonetheless, there were 5 out of 10 results that have a probability of 0.2 or greater, all of which were “q0q1q2q3q4 = 01010.” Note again that other (noisy) sequences had nearly comparable values. Next, the IonQ runs were done in two ways by the simulator-generated quantum circuits with step sizes of 5 and 1. The number of batches was 2, where the number of shots per batch was increased to 500 (see SI for the results). In the case of step size of 1, “01010” was obtained at a high count (46 in total 500; or 0.092 as probability value), corresponding to the simulator case. However, such a correspondence was not obtained in the case of step size of 5. For the IonQ system which would be susceptible to noises, care should be necessary for the fact that multiple undesirable qubit sequences appeared with sizable counts, especially in the case of deeper quantum circuits (step size of 5). Note that reported average one-qubit and two-qubit gate fidelities of the present IonQ hardware [27] are 99.5% and 97.5%, respectively, and SPAM (state preparation and measurement) error rate is 0.7%. It can thus be guessed that a boost of the number of CNOT gates due to the increase of many-body interaction terms is responsible for the poor results.

The simulation results (listed in SI) by the step size of 5 were fairly acceptable even for Trial 3, where 3 out of 5 results were the correct “01011” sequence with the highest probability. In contrast, no valid result was obtainable in the case of step size of 1 in which there were several undesired sequences with high probabilities. The IonQ runs with the quantum circuits of step size of 5 were done in two batches (each 500 shots), and no correspondence to the simulator results was found. Further IonQ runs of step size of 1 were additionally performed two times (each 1000 shots), but no improvements were observed unfortunately. In contrast to Trial 2, IonQ could not work enough for both step sizes, probably due to system noises (refer also to analyses for Trial 2).

4 SUMMARY

In this successive article of Ref [25], we described in detail the folding simulations of PSVKMA peptide on a two-dimensional lattice model [7, 8] by using the blueqat AutoQML simulator [20], by which target QUBO [17] problems were converted to QAOA [21, 22] ones automatically. The resulted quantum circuits could be passed to the IonQ quantum system [24, 27]. Three types of folding were tried. The first trial (Trial 1) of three qubits were successful by both the blueqat simulator and IonQ. However, the situations became harder for the second and third trials (Trials 2 and 3) of five qubits, especially for IonQ having a vulnerability to noises. Anyhow, a good usability of web-based quantum computation environment by the blueqat AutoQML was shown for the peptide folding problem. Note here that Sarkar et al. [28] recently reported a pioneering work by a “QUBO using QAOA” approach for DNA sequence reconstructions on the D-Wave system [16].

When hardware developments proceed from noisy intermediate-scale type quantum computer (NISQ) to fault-tolerant type quantum computer (FTQC), three-dimensional folding [10] for more realistic proteins is expectable, beside incremental costs for the actual usages. Such a future realization of quantum computation-based folding may be complementary with recent machine learning-based folding systems [29, 30]. Another route of improvements for folding may be a refinement of residue-residue interaction parameters taken from CG dynamics simulations [31, 32].

Acknowledgments

Y. M. is grateful for a support by Rikkyo SFR. K. S. acknowledges the support from the JST PRESTO project “Quantum Software” (JPMJPR1914).

REFERENCES
 
© 2023 Society of Computer Chemistry, Japan

This is an open-access article distributed under the terms of the Creative Commons Attribution Non-Commercial No Derivatives (CC BY-NC-ND) 4.0 License.
https://creativecommons.org/licenses/by-nc-nd/4.0/deed.ja
feedback
 Top