ProteinFoldingResult, XYZ files and plotting of folded Proteins for ProteinFoldingProblem.

qiskit_nature/problems/sampling/protein_folding/protein_folding_problem.py

@@ -1,6 +1,6 @@
 # This code is part of Qiskit.
 #
-# (C) Copyright IBM 2021.
+# (C) Copyright IBM 2021,2022.
 #
 # This code is licensed under the Apache License, Version 2.0. You may
 # obtain a copy of this license in the LICENSE.txt file in the root directory
@@ -10,9 +10,14 @@
 # copyright notice, and modified files need to carry a notice indicating
 # that they have been altered from the originals.
 """Defines a protein folding problem that can be passed to algorithms."""
+from __future__ import annotations
+
 from typing import Union, List
 
 from qiskit.opflow import PauliSumOp, PauliOp
+from qiskit.algorithms import MinimumEigensolverResult
+
+import qiskit_nature.results.protein_folding_result as pfr
 
 from .peptide.peptide import Peptide
 from .interactions.interaction import Interaction
@@ -88,12 +93,21 @@ def _qubit_op_full(self) -> Union[PauliOp, PauliSumOp]:
         qubit_operator = self._qubit_op_builder._build_qubit_op()
         return qubit_operator
 
-    # TODO will be implemented in another issue, including the type hint
-    def interpret(self):
-        pass
+    def interpret(self, raw_result: MinimumEigensolverResult) -> pfr.ProteinFoldingResult:
+        """
+        Returns a ProteinFoldingResult object that will allow us to interpret the result obtained.
+        For now we are only interested in the sequence with the biggest amplitude from the eigenstate.
+        """
+        best_sequence = max(raw_result.eigenstate, key=raw_result.eigenstate.get)
+        return pfr.ProteinFoldingResult(self, best_sequence)
 
     @property
     def unused_qubits(self) -> List[int]:
         """Returns the list of indices for qubits in the original problem formulation that were
         removed during compression."""
         return self._unused_qubits
+
+    @property
+    def peptide(self) -> Peptide:
+        """Returns the peptide defining the protein subject to the folding problem."""
+        return self._peptide

qiskit_nature/problems/sampling/sampling_problem.py

@@ -14,6 +14,9 @@
 from typing import Union
 
 from qiskit.opflow import PauliSumOp, PauliOp
+from qiskit.algorithms import MinimumEigensolverResult
+
+from qiskit_nature.results import EigenstateResult
 
 
 class SamplingProblem(ABC):
@@ -26,6 +29,6 @@ def qubit_op(self) -> Union[PauliOp, PauliSumOp]:
 
     # TODO type hint will be addressed later on
     @abstractmethod
-    def interpret(self):
+    def interpret(self, raw_result: MinimumEigensolverResult) -> EigenstateResult:
         """Interprets results of an optimization."""
         pass

qiskit_nature/results/protein_folding_result.py

@@ -1,6 +1,6 @@
 # This code is part of Qiskit.
 #
-# (C) Copyright IBM 2021.
+# (C) Copyright IBM 2021,2022.
 #
 # This code is licensed under the Apache License, Version 2.0. You may
 # obtain a copy of this license in the LICENSE.txt file in the root directory
@@ -10,32 +10,70 @@
 # copyright notice, and modified files need to carry a notice indicating
 # that they have been altered from the originals.
 """The protein folding result."""
+from __future__ import annotations
 from typing import Union
+import numpy as np
 
 from qiskit.opflow import PauliOp
 
-from qiskit_nature.problems.sampling.protein_folding.protein_folding_problem import (
-    ProteinFoldingProblem,
-)
+import qiskit_nature.problems.sampling.protein_folding.protein_folding_problem as pfp
+
 from qiskit_nature.results import EigenstateResult
+from .protein_folding_tools.protein_decoder import ProteinDecoder
+from .protein_folding_tools.protein_xyz import ProteinXYZ
+from .protein_folding_tools.protein_plotter import ProteinPlotter
 
 
 class ProteinFoldingResult(EigenstateResult):
-    """The protein folding result."""
+    """
+    The Protein Folding Result.
+    """
 
     def __init__(
         self,
-        protein_folding_problem: ProteinFoldingProblem,
+        protein_folding_problem: pfp.ProteinFoldingProblem,
         best_sequence: Union[str, PauliOp],
     ) -> None:
+        """
+        Args:
+            protein_folding_problem: The protein folding problem that led to the result.
+            best_sequence: The best sequence in the result eigenstate.
+
+        """
         super().__init__()
         self._protein_folding_problem = protein_folding_problem
-        self._best_sequence: str = best_sequence
+        self._best_sequence = best_sequence
+        self._unused_qubits = self._protein_folding_problem.unused_qubits
+        self._main_chain_lenght = len(
+            self._protein_folding_problem.peptide.get_main_chain.main_chain_residue_sequence
+        )
+        self._side_chain_hot_vector = (
+            self._protein_folding_problem.peptide.get_side_chain_hot_vector()
+        )
+
+        self._protein_decoder = ProteinDecoder(
+            best_sequence=self._best_sequence,
+            side_chain_hot_vector=self._side_chain_hot_vector,
+            fifth_bit=5 in self._unused_qubits[:6],
+        )
+
+        self._protein_xyz = ProteinXYZ(
+            self.protein_decoder.get_main_turns(),
+            self.protein_decoder.get_side_turns(),
+            self._protein_folding_problem.peptide,
+        )
+
+    @property
+    def protein_decoder(self) -> ProteinDecoder:
+        """Returns a ProteinDecoder.
+        This class will interpret the result bitstring and return the encoded information."""
+        return self._protein_decoder
 
     @property
-    def protein_folding_problem(self) -> ProteinFoldingProblem:
-        """Returns the protein folding problem."""
-        return self._protein_folding_problem
+    def protein_xyz(self) -> ProteinXYZ:
+        """Returns a ProteinXYZ. This class will take the encoded turns and
+        generate the position of every bead in the main and side chains."""
+        return self._protein_xyz
 
     @property
     def best_sequence(self) -> str:
@@ -48,7 +86,7 @@ def get_result_binary_vector(self) -> str:
         number of qubits in the original objective function. This method calculates the original
         version of the solution vector. Bits that can take any value without changing the
         solution are denoted by '*'."""
-        unused_qubits = self._protein_folding_problem.unused_qubits
+        unused_qubits = self._unused_qubits
         result = []
         offset = 0
         size = len(self._best_sequence)
@@ -60,3 +98,30 @@ def get_result_binary_vector(self) -> str:
             result.append(self._best_sequence[index])
 
         return "".join(result[::-1])
+
+    def get_xyz_file(self, name: str = None, output_data: bool = False) -> np.ndarray:
+        """
+        Generates an xyz_file and returns an array with the data in that file.
+        Args:
+            name: Name of the file to be generated.
+            output_data: Boolean indicating whether we want to generate the file or not.
+        Returns:
+            An array with the data to be stored in the xyz file.
+        """
+        if name is None:
+            name = str(self._protein_folding_problem.peptide.get_main_chain.main_chain_residue_sequence)
+        return self.protein_xyz.get_xyz_file(name, output_data)
+
+    def plot_folded_protein(
+        self, title: str = "Protein Structure", ticks: bool = True, grid: bool = False
+    ) -> None:
+        """
+        Plots the molecule in 3D
+        Args:
+            title: The title of the plot
+            ticks: Boolean for showing ticks in the graphic
+            grid: Boolean for showing the grid in the graphic
+
+        """
+        protein_plotter = ProteinPlotter(self)
+        protein_plotter.plot(title=title, ticks=ticks, grid=grid)

qiskit_nature/results/protein_folding_tools/__init__.py

@@ -0,0 +1,11 @@
+# This code is part of Qiskit.
+#
+# (C) Copyright IBM 2022.
+#
+# This code is licensed under the Apache License, Version 2.0. You may
+# obtain a copy of this license in the LICENSE.txt file in the root directory
+# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0.
+#
+# Any modifications or derivative works of this code must retain this
+# copyright notice, and modified files need to carry a notice indicating
+# that they have been altered from the originals.

qiskit_nature/results/protein_folding_tools/protein_decoder.py

@@ -0,0 +1,107 @@
+# This code is part of Qiskit.
+#
+# (C) Copyright IBM 2022.
+#
+# This code is licensed under the Apache License, Version 2.0. You may
+# obtain a copy of this license in the LICENSE.txt file in the root directory
+# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0.
+#
+# Any modifications or derivative works of this code must retain this
+# copyright notice, and modified files need to carry a notice indicating
+# that they have been altered from the originals.
+"""An auxiliary class that gets the turns in the main and side chain of a molecule
+ in ProteinFoldingResult """
+from typing import Union, List, Tuple
+
+
+class ProteinDecoder:
+
+    """This class handles the decoding of the compact solution in ProteinFoldingProblem
+    and returns the information encoded in the result about the turns
+    associated to the main and side chains.
+    """
+
+    def __init__(
+        self, best_sequence: str, side_chain_hot_vector: List[bool], fifth_bit: bool
+    ) -> None:
+        """
+        Args:
+            best_sequence: Will be the sequence decoded.
+            side_chain_hot_vector: boolean list with the position of the side chains
+            unused_qubits: list of qubits compressed in the ProteinFoldingProblem.
+        """
+        self._best_sequence = best_sequence
+        self._side_chain_hot_vector = side_chain_hot_vector
+        self._fifth_bit = fifth_bit
+        self._main_chain_lenght = len(side_chain_hot_vector)
+
+    def _bitstring_to_turns(self, bitstring: str) -> List[int]:
+        """
+        Takes a bitstring encoding the turns of a chain and retuns the turns as a list of int.
+
+        Args:
+            bitstring: string containing the encoded information.
+        Returns:
+            A list of integers decoding the bitstring.
+        """
+        bitstring = bitstring[::-1]
+        encoding = {"00": 0, "01": 1, "10": 2, "11": 3}
+        lenght_turns = len(bitstring) // 2
+        return [encoding[bitstring[2 * i : 2 * (i + 1)]] for i in range(lenght_turns)]
+
+    def _split_bitstring(self) -> Tuple[int, int]:
+        """Returns the amount of bits in the compact solution corresponding
+        to each property they encode."""
+        n_qbits_encoding_main_turns = 2 * (self._main_chain_lenght - 3) - (self._fifth_bit)
+        n_qbits_encoding_side_turns = 2 * sum(self._side_chain_hot_vector)
+        return n_qbits_encoding_main_turns, n_qbits_encoding_side_turns
+
+    def get_main_turns(self) -> List[int]:
+        """
+        Returns the list of turns for the molecule corresponding to best_sequence.
+        The first element of the list corresponds to the turn of the second aminoacid in the peptide.
+        Returns:
+                A list of integers representing the sequence of turns on the molecule
+
+        Notes:
+                The bitstring will end in 0010 corresponding to turn1=(01) and turn2=(00)
+                If the second bead doesn't have a side chain the 6th bit
+                can be set to 1 without loss of generality.
+                In that case index (5) will belong to the list of unused qubits.
+                The amount of qubits needed to encode the turns will be 2(N-3) - 1
+                if no side chain on second main bead or 2(N-3) otherwise.
+        """
+
+        main_turns_bitstring = self._best_sequence[-self._split_bitstring()[0] :] + "0010"
+
+        if self._fifth_bit:
+            main_turns_bitstring = main_turns_bitstring[:-5] + "1" + main_turns_bitstring[-5:]
+
+        return self._bitstring_to_turns(main_turns_bitstring)
+
+    def get_side_turns(self) -> List[Union[None, int]]:
+        """
+        Returns the list of turns from the main bead corresponding to the side chains.
+        None corresponds to no side chain from that main bead.
+
+        Returns:
+            A list with either an number associated to a turn from
+            the main bead or None if there is no side bead.
+
+        """
+        n, m = self._split_bitstring()
+
+        side_turns_bitstring = self._best_sequence[-n - m : -n]
+
+        side_turns = self._bitstring_to_turns(side_turns_bitstring)
+
+        result = []
+        counter = 0
+        for element in self._side_chain_hot_vector:
+            if element:
+                result.append(side_turns[counter])
+                counter += 1
+            else:
+                result.append(None)
+
+        return result