# Membrane Environment Modulates Ligand-Binding Propensity of P2Y12 Receptor

Fatemeh Haghighi, Semen Yesylevskyy, Siamak Davani, Christophe Ramseyer

► To cite this version:

Fatemeh Haghighi, Semen Yesylevskyy, Siamak Davani, Christophe Ramseyer. Membrane Environment Modulates Ligand-Binding Propensity of P2Y12 Receptor. Pharmaceutics, 2021, 13 (4), pp.524. 10.3390/pharmaceutics13040524. [https://doi.org/10.3390/pharmaceutics13040524](https://doi.org/10.3390/pharmaceutics13040524)

HAL Id: hal-04226833

https://hal.science/hal-04226833

Submitted on 28 Jun 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Membrane Environment Modulates Ligand-Binding Propensity of P2Y12 Receptor

Fatemeh Haghighi, Semen Yesylevskyy, Siamak Davani, Christophe Ramseyer

► To cite this version:

Fatemeh Haghighi, Semen Yesylevskyy, Siamak Davani, Christophe Ramseyer. Membrane Environment Modulates Ligand-Binding Propensity of P2Y12 Receptor. Pharmaceutics, 2021, 13 (4), pp.524. 10.3390/pharmaceutics13040524. [https://doi.org/10.3390/pharmaceutics13040524](https://doi.org/10.3390/pharmaceutics13040524). [hal-04226833](https:// HAL Id: )

Submitted on 28 Jun 2024

**HAL** is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire **HAL**, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_1.png" /></figure>
Article

# Membrane Environment Modulates Ligand-Binding Propensity of P2Y12 Receptor

Fatemeh Haghighi <fnref n="1" />,<fnref n="2" /> ID, Semen Yesylevskyy <fnref n="1" />,<fnref n="3" />, Siamak Davani <fnref n="2" />,<fnref n="4" /> ID and Christophe Ramseyer <fnref n="1" />,*

<figure><img src="image_2.png" /></figure>

check for
updates

**Citation:** Haghighi, F.; Yesylevskyy, S.; Davani, S.; Ramseyer, C.

Membrane Environment Modulates Ligand-Binding Propensity of P2Y12 Receptor. *Pharmaceutics* **2021**, *13*, 524.
https://doi.org/10.3390/pharmaceutics13040524

Academic Editor: Fabrizio Dal Piaz

Received: 9 March 2021

Accepted: 6 April 2021

Published: 9 April 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_3.png" /></figure>

**Copyright:** © 2021 by the authors.
Licensee MDPI, Basel, Switzerland.
This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

**Abstract:** The binding of natural ligands and synthetic drugs to the P2Y12 receptor is of great interest because of its crucial role in platelets activation and the therapy of arterial thrombosis. Up to now, all computational studies of P2Y12 concentrated on the available crystal structures, while the role of intrinsic protein dynamics and the membrane environment in the functioning of P2Y12 was not clear. In this work, we performed all-atom molecular dynamics simulations of the full-length P2Y₁₂ receptor in three different membrane environments and in two possible conformations derived from available crystal structures. The binding of ticagrelor, its two major metabolites, adenosine diphosphate (ADP) and 2-Methylthioadenosine diphosphate (2MeS-ADP) as agonist, and ethyl 6-[4-(benzylsulfonylcarbamoyl)piperidin-1-yl]-5-cyano-2-methylpyridine-3-carboxylate (AZD1283) as antagonist were assessed systematically by means of ensemble docking. It is shown that the binding of all ligands becomes systematically stronger with the increase of the membrane rigidity. Binding of all ligands to the agonist-bound-like conformations is systematically stronger in comparison to antagonist-bound-likes ones. This is dramatically opposite to the results obtained for static crystal structures. Our results show that accounting for internal protein dynamics, strongly modulated by its lipid environment, is crucial for correct assessment of the ligand binding to P2Y12.

**Keywords:** ticagrelor; P2Y₁₂ receptors; platelets; molecular dynamics

## 1. Introduction

Platelets are major drug targets in case of thromboembolic disorders. Pharmaceutical inhibition of platelet aggregation is necessary for preventing complications after an acute coronary syndrome (ACS) and the formation and progression of thrombotic processes. Adenosine diphosphate (ADP) is an important signal molecule for platelet function, which is involved in their physiological and pathological responses [1]. The inhibition of ADP-induced platelet aggregation is an effective approach for treating thrombotic events in clinical practice. ADP binds to purinergic receptors P2Y1 and P2Y12 which are seven-transmembrane (7TM) G-protein-coupled receptor (GPCR) located in the platelet cell membrane [2]. These receptors can be activated by extracellular nucleotides and trigger the completion of the platelet aggregation [3,4]. Therefore, ADP receptors are nowadays considered as the main targets for antiplatelet agents [5]. The P2Y12 receptor is abundantly expressed in the platelet membrane [6] and is a key player in primary hemostasis and arterial thrombosis, which makes it one of the most prominent drug targets for the inhibition of platelet aggregation [5,6].

<footnote n="1">¹ Laboratoire Chrono Environnement UMR CNRS 6249, Université de Bourgogne Franche-Comté, 16 route de Gray, CEDEX, 25030 Besançon, France; fateme.haghighi@hotmail.com (F.H.); yesint4@gmail.com (S.Y.)</footnote>

<footnote n="2">² Unité de Recherche EA 3920, Université de Bourgogne Franche-Comté, 19 rue Ambroise Paré, CEDEX, 25000 Besançon, France; siamak.davani@univ-fcomte.fr</footnote>

<footnote n="3">³ Department of Physics of Biological Systems, Institute of Physics of the National Academy of Sciences of Ukraine, Prospect Nauky 46, 03028 Kyiv, Ukraine</footnote>

<footnote n="4">⁴ Laboratoire de Pharmacologie Clinique et Toxicologie-CHU de Besançon, 3 boulevard A. Fleming, CEDEX, 25030 Besançon, France</footnote>

<footnote>* Correspondence: christophe.ramseyer@univ-fcomte.fr; Tel.: +33-381666482</footnote>
Article

# Membrane Environment Modulates Ligand-Binding Propensity of P2Y12 Receptor

Fatemeh Haghighi <fnref n="1" />,<fnref n="2" /> ID, Semen Yesylevskyy <fnref n="1" />,<fnref n="3" />, Siamak Davani <fnref n="2" />,<fnref n="4" /> ID and Christophe Ramseyer <fnref n="1" />,*

<figure><img src="image_2.png" /></figure>

check for
updates

**Citation:** Haghighi, F.; Yesylevskyy, S.; Davani, S.; Ramseyer, C.

Membrane Environment Modulates Ligand-Binding Propensity of P2Y12 Receptor. *Pharmaceutics* **2021**, *13*, 524.
https://doi.org/10.3390/pharmaceutics13040524

Academic Editor: Fabrizio Dal Piaz

Received: 9 March 2021

Accepted: 6 April 2021

Published: 9 April 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_3.png" /></figure>

**Copyright:** © 2021 by the authors.
Licensee MDPI, Basel, Switzerland.
This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

**Abstract:** The binding of natural ligands and synthetic drugs to the P2Y12 receptor is of great interest because of its crucial role in platelets activation and the therapy of arterial thrombosis. Up to now, all computational studies of P2Y12 concentrated on the available crystal structures, while the role of intrinsic protein dynamics and the membrane environment in the functioning of P2Y12 was not clear. In this work, we performed all-atom molecular dynamics simulations of the full-length P2Y₁₂ receptor in three different membrane environments and in two possible conformations derived from available crystal structures. The binding of ticagrelor, its two major metabolites, adenosine diphosphate (ADP) and 2-Methylthioadenosine diphosphate (2MeS-ADP) as agonist, and ethyl 6-[4-(benzylsulfonylcarbamoyl)piperidin-1-yl]-5-cyano-2-methylpyridine-3-carboxylate (AZD1283) as antagonist were assessed systematically by means of ensemble docking. It is shown that the binding of all ligands becomes systematically stronger with the increase of the membrane rigidity. Binding of all ligands to the agonist-bound-like conformations is systematically stronger in comparison to antagonist-bound-likes ones. This is dramatically opposite to the results obtained for static crystal structures. Our results show that accounting for internal protein dynamics, strongly modulated by its lipid environment, is crucial for correct assessment of the ligand binding to P2Y12.

**Keywords:** ticagrelor; P2Y₁₂ receptors; platelets; molecular dynamics

## 1. Introduction

Platelets are major drug targets in case of thromboembolic disorders. Pharmaceutical inhibition of platelet aggregation is necessary for preventing complications after an acute coronary syndrome (ACS) and the formation and progression of thrombotic processes. Adenosine diphosphate (ADP) is an important signal molecule for platelet function, which is involved in their physiological and pathological responses [1]. The inhibition of ADP-induced platelet aggregation is an effective approach for treating thrombotic events in clinical practice. ADP binds to purinergic receptors P2Y1 and P2Y12 which are seven-transmembrane (7TM) G-protein-coupled receptor (GPCR) located in the platelet cell membrane [2]. These receptors can be activated by extracellular nucleotides and trigger the completion of the platelet aggregation [3,4]. Therefore, ADP receptors are nowadays considered as the main targets for antiplatelet agents [5]. The P2Y12 receptor is abundantly expressed in the platelet membrane [6] and is a key player in primary hemostasis and arterial thrombosis, which makes it one of the most prominent drug targets for the inhibition of platelet aggregation [5,6].

<footnote n="1">¹ Laboratoire Chrono Environnement UMR CNRS 6249, Université de Bourgogne Franche-Comté, 16 route de Gray, CEDEX, 25030 Besançon, France; fateme.haghighi@hotmail.com (F.H.); yesint4@gmail.com (S.Y.)</footnote>

<footnote n="2">² Unité de Recherche EA 3920, Université de Bourgogne Franche-Comté, 19 rue Ambroise Paré, CEDEX, 25000 Besançon, France; siamak.davani@univ-fcomte.fr</footnote>

<footnote n="3">³ Department of Physics of Biological Systems, Institute of Physics of the National Academy of Sciences of Ukraine, Prospect Nauky 46, 03028 Kyiv, Ukraine</footnote>

<footnote n="4">⁴ Laboratoire de Pharmacologie Clinique et Toxicologie-CHU de Besançon, 3 boulevard A. Fleming, CEDEX, 25030 Besançon, France</footnote>

<footnote>* Correspondence: christophe.ramseyer@univ-fcomte.fr; Tel.: +33-381666482</footnote>
Currently, two classes of P2Y12 antagonists have been developed and approved by the FDA as antiplatelet agents, namely the irreversible thienopyridines and the reversible ATP analogues. Thienopyridine compounds such as ticlopidine, clopidogrel, and prasugrel irreversibly inhibit P2Y12 receptor [7]. They are prodrugs that must be metabolized in the liver [8]. The active metabolite of clopidogrel and prasugrel covalently binds to cysteine residues of P2Y12 precluding the binding of ADP [9,10]. Reversible P2Y12 receptor antagonists are modeled after ATP as a scaffold [11]. Ticagrelor is the compound from this group that is widely used in clinical practice [11,12]. It binds to the receptor directly without the need for hepatic bioactivation [13], which results in a faster onset of action, less interindividual variation in effectiveness, and greater control of platelet aggregation inhibition [14]. It is reported that the primary metabolite of ticagrelor is also active with potency equal to the parent compound [12].

Limitations of existing drugs, such as a very long half-life of clopidogrel and adverse effects of ticagrelor [11,15], indicate that there is an unfulfilled demand for a new generation of P2Y12 inhibitors. Unfortunately, ligand recognition and the function of P2Y12 remain poorly understood at the molecular level due to unclear or even conflicting structure-activity relationship (SAR) data [16,17], which hampers the development of new inhibitors. It is well known that the lipids regulate the location and activity of many membrane proteins as well as forming lipid microdomains that control the localization and interactions of proteins involved in cell signaling [18,19]. Being an integral membrane protein, P2Y12 receptor is also likely to be influenced and regulated by its lipid environment. Particularly, the impact of ticagrelor on the lipid composition of platelet plasma membrane [20] and the role of cholesterol on the activation of P2Y12 receptor was recently studied [21]. It was shown that the location of P2Y12 in ordered lipid microdomains ("rafts") is required during ADP-mediated platelet activation [21]. Moreover, the active metabolite of clopidogrel breaks up the homooligomers of P2Y12 receptor into non-functional dimers and monomers that are located outside of lipid rafts [10]. The effect of ADP and ticagrelor on the detergent-resistant membrane (DRM) was also examined recently [22].

Experimental studies of the influence of lipid environment on membrane proteins are complex, tedious, and time consuming, which recently boosted the development of computational approaches to this problem. In the case of P2Y12, early numerical studies [23,24] were pushed forward by resolution of the crystallographic structures of the P2Y12 complexed with the full agonist 2-methylthioadenosine-5'-diphosphate (2MeSADP, PDB ID: 4PXZ), with the partial agonist 2-methylthio-adenosine-5'-triphosphate (2MeSATP, PDB ID: 4PY0) [25] and with the antagonist ethyl 6-(4-((benzylsulfonyl)carbamoyl)piperidine-1-yl)-5-cyano-2-methylnicotinate (AZD1283, PDB ID: 4NTJ) [25]. These structures differ substantially by the conformation of transmembrane helix 6 (TM6). In the agonist-bound 4PXZ, this helix is bent inwards toward the active site (this structure will be referred as "closed" hereafter). In antagonist-bound 4NTJ, this helix is straight, and the active site is significantly more exposed (this structure will be called "open" hereafter). These large-scale rearrangements in the highly malleable extracellular regions are required to allow agonist accessing the binding pocket. The previous structural and docking studies [26] were not able to provide valuable insight into the pharmacology and mechanisms of action of agonists and different classes of antagonists because they were performed before these structural differences become known. For instance, ticagrelor could not be docked into a pose similar to one of 2MeS-ADP because the presence of bulky N6 substituents would cause a steric clash with the helixes V and VI [26].

Recently, an extensive molecular docking study of the major classes of substances, previously reported as P2Y12 ligands, was performed [27]. The authors tried to rationalize the main SAR findings previously reported for each class of the ligands. It was shown that antagonists such as ureas, sulfonylureas, sulfonamides, anthraquinones, and glutamic acid piperazines docked readily to the antagonist-bound P2Y12 structure, while various nucleotide derivatives docked readily to the agonist-bound structure. However, they were unable to dock ticagrelor to the agonist-bound P2Y12 structure. The "hybrid" receptor
Currently, two classes of P2Y12 antagonists have been developed and approved by the FDA as antiplatelet agents, namely the irreversible thienopyridines and the reversible ATP analogues. Thienopyridine compounds such as ticlopidine, clopidogrel, and prasugrel irreversibly inhibit P2Y12 receptor [7]. They are prodrugs that must be metabolized in the liver [8]. The active metabolite of clopidogrel and prasugrel covalently binds to cysteine residues of P2Y12 precluding the binding of ADP [9,10]. Reversible P2Y12 receptor antagonists are modeled after ATP as a scaffold [11]. Ticagrelor is the compound from this group that is widely used in clinical practice [11,12]. It binds to the receptor directly without the need for hepatic bioactivation [13], which results in a faster onset of action, less interindividual variation in effectiveness, and greater control of platelet aggregation inhibition [14]. It is reported that the primary metabolite of ticagrelor is also active with potency equal to the parent compound [12].

Limitations of existing drugs, such as a very long half-life of clopidogrel and adverse effects of ticagrelor [11,15], indicate that there is an unfulfilled demand for a new generation of P2Y12 inhibitors. Unfortunately, ligand recognition and the function of P2Y12 remain poorly understood at the molecular level due to unclear or even conflicting structure-activity relationship (SAR) data [16,17], which hampers the development of new inhibitors. It is well known that the lipids regulate the location and activity of many membrane proteins as well as forming lipid microdomains that control the localization and interactions of proteins involved in cell signaling [18,19]. Being an integral membrane protein, P2Y12 receptor is also likely to be influenced and regulated by its lipid environment. Particularly, the impact of ticagrelor on the lipid composition of platelet plasma membrane [20] and the role of cholesterol on the activation of P2Y12 receptor was recently studied [21]. It was shown that the location of P2Y12 in ordered lipid microdomains ("rafts") is required during ADP-mediated platelet activation [21]. Moreover, the active metabolite of clopidogrel breaks up the homooligomers of P2Y12 receptor into non-functional dimers and monomers that are located outside of lipid rafts [10]. The effect of ADP and ticagrelor on the detergent-resistant membrane (DRM) was also examined recently [22].

Experimental studies of the influence of lipid environment on membrane proteins are complex, tedious, and time consuming, which recently boosted the development of computational approaches to this problem. In the case of P2Y12, early numerical studies [23,24] were pushed forward by resolution of the crystallographic structures of the P2Y12 complexed with the full agonist 2-methylthioadenosine-5'-diphosphate (2MeSADP, PDB ID: 4PXZ), with the partial agonist 2-methylthio-adenosine-5'-triphosphate (2MeSATP, PDB ID: 4PY0) [25] and with the antagonist ethyl 6-(4-((benzylsulfonyl)carbamoyl)piperidine-1-yl)-5-cyano-2-methylnicotinate (AZD1283, PDB ID: 4NTJ) [25]. These structures differ substantially by the conformation of transmembrane helix 6 (TM6). In the agonist-bound 4PXZ, this helix is bent inwards toward the active site (this structure will be referred as "closed" hereafter). In antagonist-bound 4NTJ, this helix is straight, and the active site is significantly more exposed (this structure will be called "open" hereafter). These large-scale rearrangements in the highly malleable extracellular regions are required to allow agonist accessing the binding pocket. The previous structural and docking studies [26] were not able to provide valuable insight into the pharmacology and mechanisms of action of agonists and different classes of antagonists because they were performed before these structural differences become known. For instance, ticagrelor could not be docked into a pose similar to one of 2MeS-ADP because the presence of bulky N6 substituents would cause a steric clash with the helixes V and VI [26].

Recently, an extensive molecular docking study of the major classes of substances, previously reported as P2Y12 ligands, was performed [27]. The authors tried to rationalize the main SAR findings previously reported for each class of the ligands. It was shown that antagonists such as ureas, sulfonylureas, sulfonamides, anthraquinones, and glutamic acid piperazines docked readily to the antagonist-bound P2Y12 structure, while various nucleotide derivatives docked readily to the agonist-bound structure. However, they were unable to dock ticagrelor to the agonist-bound P2Y12 structure. The "hybrid" receptor
resembling the agonist-bound P2Y12 except for the top portion of TM6, which was taken
from the antagonist-bound P2Y12 structure, was necessary for successful ticagrelor docking.
This example shows the limitations of available P2Y12 models and emphasizes a need
for more elaborate computational studies, which take into account protein mobility and
conformational changes in its native membrane environment.

In the present study, we propose the most realistic computational model of P2Y12
available to date. First, we consider the receptor in its native lipid bilayer environment. The
dynamics of the receptor is investigated as a function of bilayer composition. Phosphatidyl-
choline (PC), sphingomyelin (SM), and SM/cholesterol mixture (referred as "raft" hereafter)
environments are studied. Second, both agonist and antagonist-bound crystal structures
are used as initial conformations in each lipid environment. Third, an extensive ensemble
docking analysis is made on a representative set of simulation snapshots extracted from
equilibrated molecular dynamics (MD) trajectories, which takes into account both the local
flexibility of the binding site and the large-scale protein motions. This methodology was
already used with great success for the proteins with complex internal dynamics, such as
STAT1/STAT3 [28,29] and human serum albumin [30,31].

We demonstrate that the lipid environment imposes significant constraints into the
dynamics and internal flexibility of P2Y12, which in turn influences its binding to a variety
of ligands. In general, the more ordered the lipid environment, the higher the binding
propensity of P2Y12 to all studied ligands. The bulky nucleotide ligands, such as ticagrelor,
can access the binding pockets of both agonist and antagonist-bound P2Y12 conformations
with the help of internal protein dynamics.

**2. Materials and Methods**

**2.1. Homology Modeling**

In this work, we considered two crystal structures of P2Y12: 4PXZ:A in complex with 2MeS-ADP agonist [26] (referred as “closed” conformation) and 4NTJ:A in complex with AZD1283 antagonist drug [25] (referred as “open” conformation).

The closed structure of P2Y12 (PDB code 4PXZ:A) was chosen as the starting point for
homology modeling. It covers only the residues 15–305, which constitute the major part of
transmembrane core of this membrane protein. The crystallized structure is an engineered
chimeric protein with soluble cytochrome B562 inserted between the residues 223 and 224
of P2Y12. According to UniProt annotation (UniProtKB: Q9H244), the secondary structure
of the missing regions 1–14 and 306–342 is not known and presumed to be unstructured.

The homology modeling approach was used to reconstruct the full-length structure of P2Y12. The cytochrome insertion was removed, and the GalaxyTBM server [32] was used for homology modeling based on the full-length sequence. There are no suitable homology templates for the long C-terminal region 306–342; thus, the structure generated for this region was not satisfactory. Particularly, it protruded deep into the membrane, while the C-terminus should be located in the cytoplasm according to UniProt annotation. That is why we performed de novo structure prediction for the region 302–342 using the QUARK server [33]. Then, the predicted structure was aligned with the rest of the protein using the overlapping helical region 302–305 and manually attached to the transmembrane segment.

In order to construct the open structure of the protein, we first aligned the obtained
homology model of the closed form with the crystal structure 4NTJ:A and then transferred
all atomic coordinates from the crystal structure to our model. This resulted in the models
of closed and open conformations which share the same reconstructed peripheral regions
missed in the crystal structures.

**2.2. *Binding Pocket Identification***

According to the closed crystal structure, the residues reported to bind with 2MES-ADP are ARG19, ARG93, CYS97, SER101, VAL102, TYR105, PHE106, TYR109, MET152, LEU155, SER156, ASN159, THR163, CYS175, LYS179, HIS187, VAL190, ASN191, CYS194, ARG256, TYR259, GLN263, and LYS280 [26]. In the open crystal structure, most of these
resembling the agonist-bound P2Y12 except for the top portion of TM6, which was taken
from the antagonist-bound P2Y12 structure, was necessary for successful ticagrelor docking.
This example shows the limitations of available P2Y12 models and emphasizes a need
for more elaborate computational studies, which take into account protein mobility and
conformational changes in its native membrane environment.

In the present study, we propose the most realistic computational model of P2Y12
available to date. First, we consider the receptor in its native lipid bilayer environment. The
dynamics of the receptor is investigated as a function of bilayer composition. Phosphatidyl-
choline (PC), sphingomyelin (SM), and SM/cholesterol mixture (referred as "raft" hereafter)
environments are studied. Second, both agonist and antagonist-bound crystal structures
are used as initial conformations in each lipid environment. Third, an extensive ensemble
docking analysis is made on a representative set of simulation snapshots extracted from
equilibrated molecular dynamics (MD) trajectories, which takes into account both the local
flexibility of the binding site and the large-scale protein motions. This methodology was
already used with great success for the proteins with complex internal dynamics, such as
STAT1/STAT3 [28,29] and human serum albumin [30,31].

We demonstrate that the lipid environment imposes significant constraints into the
dynamics and internal flexibility of P2Y12, which in turn influences its binding to a variety
of ligands. In general, the more ordered the lipid environment, the higher the binding
propensity of P2Y12 to all studied ligands. The bulky nucleotide ligands, such as ticagrelor,
can access the binding pockets of both agonist and antagonist-bound P2Y12 conformations
with the help of internal protein dynamics.

**2. Materials and Methods**

**2.1. Homology Modeling**

In this work, we considered two crystal structures of P2Y12: 4PXZ:A in complex with 2MeS-ADP agonist [26] (referred as “closed” conformation) and 4NTJ:A in complex with AZD1283 antagonist drug [25] (referred as “open” conformation).

The closed structure of P2Y12 (PDB code 4PXZ:A) was chosen as the starting point for
homology modeling. It covers only the residues 15–305, which constitute the major part of
transmembrane core of this membrane protein. The crystallized structure is an engineered
chimeric protein with soluble cytochrome B562 inserted between the residues 223 and 224
of P2Y12. According to UniProt annotation (UniProtKB: Q9H244), the secondary structure
of the missing regions 1–14 and 306–342 is not known and presumed to be unstructured.

The homology modeling approach was used to reconstruct the full-length structure of P2Y12. The cytochrome insertion was removed, and the GalaxyTBM server [32] was used for homology modeling based on the full-length sequence. There are no suitable homology templates for the long C-terminal region 306–342; thus, the structure generated for this region was not satisfactory. Particularly, it protruded deep into the membrane, while the C-terminus should be located in the cytoplasm according to UniProt annotation. That is why we performed de novo structure prediction for the region 302–342 using the QUARK server [33]. Then, the predicted structure was aligned with the rest of the protein using the overlapping helical region 302–305 and manually attached to the transmembrane segment.

In order to construct the open structure of the protein, we first aligned the obtained
homology model of the closed form with the crystal structure 4NTJ:A and then transferred
all atomic coordinates from the crystal structure to our model. This resulted in the models
of closed and open conformations which share the same reconstructed peripheral regions
missed in the crystal structures.

**2.2. *Binding Pocket Identification***

According to the closed crystal structure, the residues reported to bind with 2MES-ADP are ARG19, ARG93, CYS97, SER101, VAL102, TYR105, PHE106, TYR109, MET152, LEU155, SER156, ASN159, THR163, CYS175, LYS179, HIS187, VAL190, ASN191, CYS194, ARG256, TYR259, GLN263, and LYS280 [26]. In the open crystal structure, most of these
residues (SER101, VAL102, TYR105, PHE106, TYR109, MET152, LEU155, SER156, ASN159, HIS187, VAL190, ASN191, CYS194, ARG256, TYR259, LYS280) also bind to AZD1283 [25]. Some other residues such as GLN195, PHE252, ALA255, THR260, LEU276, VAL279, and THR283 were also involved in binding [25]. The center of masses of all these residues was used in this work to determine the center of the binding pocket.

## 2.3. *The Ligands*

We studied ticagrelor (TIC) and its two major metabolites [34]: an active metabolite AR-C124910XX circulating in the blood (referred as “M8” following the naming in [12]) and the metabolite AR-C133913XX found in the urine (referred as “M5” following the naming in [12]). We also studied ADP, 2MeSADP, and the AZD1283 for comparison with ticagrelor and its metabolites.

## 2.4. *Molecular Dynamics Simulations*

In this work, the open and closed forms of P2Y12 were inserted into three different lipid environments, namely PC (palmitoyl-oleoyl-phosphocholine, POPC 16:0/18:1, SM (PSM, d18:1/16:0 palmitoylsphingomyelin, and 1:1 SM/cholesterol mixture, which mimics cholesterol-enriched microdomains also widely known as “rafts”. Note that the lengths of fatty acid chains (d18:1/16:0) were reported to be the most abundant ones in platelets membranes [22]. The membranes with inserted protein were generated with the CHARMM GUI membrane builder tool [35,36]. Na⁺ and Cl⁻ ions were added corresponding to ionic strength of ≈0.15 mol/L and adjusted to counterbalance the net charge of the protein. All systems were hydrated to approximately 50 water molecules per lipid. The content of the resulting systems is shown in Table 1.

<table><caption>**Table 1.** Lipid and ion content of the studied systems.</caption><thead><tr><td>System</td><td>Lipids</td><td>Na<sup>+</sup></td><td>Cl<sup>-</sup></td></tr></thead><tbody><tr><td>PC</td><td>POPC, 220</td><td>20</td><td>39</td></tr><tr><td>SM</td><td>PSM, 220</td><td>18</td><td>37</td></tr><tr><td>Raft</td><td>PSM, 112; Cholesterol, 112</td><td>21</td><td>40</td></tr></tbody></table>

The membrane ordering and fluidity was assessed by computing lateral diffusion coefficients of the lipids in all studied systems. The results are shown in Table S1 (see Supplementary Materials). It is clearly seen that the PC membrane is an order of magnitude more fluid than SM and raft membranes. In turn, the raft membrane is only slightly less fluid than the SM one. The diffusion coefficients were computed using the mean square displacement (RMS) method as implemented in Gromacs [37] version 5.1.2 by manual fitting of the linear part of the RMS curve.

The open form of the protein was built by substituting the transmembrane core of the closed form with the coordinates from 4NTJ:A crystal structure as described above. The open forms were built separately for three lipid environments from the corresponding pre-equilibrated closed forms.

For all systems, the following equilibration protocol was applied. First, all protein atoms were restrained, and the lipids and solvent were equilibrated for 100 ns. Then, the protein backbone was restrained, and the side chains were equilibrated with the rest of the system for another 100 ns. Finally, all restraints were removed, and the production run of at least 200 ns was performed. The equilibration was monitored by the Root Mean Square Deviation (RMSD) of the transmembrane part of the protein (residues 20–300) excluding the termini, which were fluctuating significantly during the simulations. Simulations that did not converge in 200 ns were prolonged until the RMSD stabilizes.

All MD simulations were performed in Gromacs [37] version 5.1.2 in NPT ensemble at the pressure of 1 atm maintained by Parrinello-Rahman barostat [38] with semi-isotropic pressure coupling. The Verlet cutoff scheme was used [39]. Force-switch cut-off of the
residues (SER101, VAL102, TYR105, PHE106, TYR109, MET152, LEU155, SER156, ASN159, HIS187, VAL190, ASN191, CYS194, ARG256, TYR259, LYS280) also bind to AZD1283 [25]. Some other residues such as GLN195, PHE252, ALA255, THR260, LEU276, VAL279, and THR283 were also involved in binding [25]. The center of masses of all these residues was used in this work to determine the center of the binding pocket.

## 2.3. The Ligands

We studied ticagrelor (TIC) and its two major metabolites [34]: an active metabolite AR-C124910XX circulating in the blood (referred as "M8" following the naming in [12]) and the metabolite AR-C133913XX found in the urine (referred as "M5" following the naming in [12]). We also studied ADP, 2MeSADP, and the AZD1283 for comparison with ticagrelor and its metabolites.

## 2.4. Molecular Dynamics Simulations

In this work, the open and closed forms of P2Y12 were inserted into three different lipid environments, namely PC (palmitoyl-oleoyl-phosphocholine, POPC 16:0/18:1, SM (PSM, d18:1/16:0 palmitoylsphingomyelin, and 1:1 SM/cholesterol mixture, which mimics cholesterol-enriched microdomains also widely known as "rafts". Note that the lengths of fatty acid chains (d18:1/16:0) were reported to be the most abundant ones in platelets membranes [22]. The membranes with inserted protein were generated with the CHARMM GUI membrane builder tool [35,36]. Na⁺ and Cl⁻ ions were added corresponding to ionic strength of ≈0.15 mol/L and adjusted to counterbalance the net charge of the protein. All systems were hydrated to approximately 50 water molecules per lipid. The content of the resulting systems is shown in Table 1.

<table><caption>**Table 1.** Lipid and ion content of the studied systems.</caption><thead><tr><td>System</td><td>Lipids</td><td>Na<sup>+</sup></td><td>Cl<sup>-</sup></td></tr></thead><tbody><tr><td>PC</td><td>POPC, 220</td><td>20</td><td>39</td></tr><tr><td>SM</td><td>PSM, 220</td><td>18</td><td>37</td></tr><tr><td>Raft</td><td>PSM, 112; Cholesterol, 112</td><td>21</td><td>40</td></tr></tbody></table>

The membrane ordering and fluidity was assessed by computing lateral diffusion coefficients of the lipids in all studied systems. The results are shown in Table S1 (see Supplementary Materials). It is clearly seen that the PC membrane is an order of magnitude more fluid than SM and raft membranes. In turn, the raft membrane is only slightly less fluid than the SM one. The diffusion coefficients were computed using the mean square displacement (RMS) method as implemented in Gromacs [37] version 5.1.2 by manual fitting of the linear part of the RMS curve.

The open form of the protein was built by substituting the transmembrane core of the closed form with the coordinates from 4NTJ:A crystal structure as described above. The open forms were built separately for three lipid environments from the corresponding pre-equilibrated closed forms.

For all systems, the following equilibration protocol was applied. First, all protein atoms were restrained, and the lipids and solvent were equilibrated for 100 ns. Then, the protein backbone was restrained, and the side chains were equilibrated with the rest of the system for another 100 ns. Finally, all restraints were removed, and the production run of at least 200 ns was performed. The equilibration was monitored by the Root Mean Square Deviation (RMSD) of the transmembrane part of the protein (residues 20–300) excluding the termini, which were fluctuating significantly during the simulations. Simulations that did not converge in 200 ns were prolonged until the RMSD stabilizes.

All MD simulations were performed in Gromacs [37] version 5.1.2 in NPT ensemble at the pressure of 1 atm maintained by Parrinello-Rahman barostat [38] with semi-isotropic pressure coupling. The Verlet cutoff scheme was used [39]. Force-switch cut-off of the
Van der Waals interactions was used in the region between 1.0 and 1.2 nm. Long-range electrostatics was computed with the Particle-Mesh-Ewald (PME) method [40] with the cut-off of explicit short-range electrostatic interactions at 1.2 nm. Velocity rescale thermostat [41] was used at the temperature of 310K. The CHARMM36 force field [42] was used for all components of the system. An integration step of 2 fs was used in all simulations with the bonds to hydrogen atoms converted to rigid constraints. Analysis was performed by the custom scripts based on the Pteros molecular modeling library [43].

## 2.5. Ensemble Docking Simulations

In order to account for protein conformational flexibility, multiple conformations from MD trajectories were extracted and used in docking simulations. The last 50 ns of MD trajectories were used to extract 834 equally spaced frames. For each frame, the docking of all ligands was performed. The docking volume was centered in the center of mass of the binding site and had dimensions 2 × 2 × 2 nm. Quick Vina 2 [44] docking software was used with default scoring function. The best docking pose was kept in each simulation. The methodology and scripts for ensemble docking used in our previous works [29–31] were used. No side chains were treated as flexible, because the possible protein flexibility was already accounted for by means of MD simulations.

# 3. Results and Discussion

## 3.1. Protein Structure and Flexibility

Both open and closed forms of the protein were found to be stable in all membrane environments (Figure 1). However, some differences are worth to be noted, as evident from the analysis of RMSD plots (Supplementary Figure S1).

<figure><img src="image_4.png" /><figcaption>**Figure 1.** Snapshots of equilibrated closed form of protein in phosphatidylcholine (PC) (A), sphingomyelin (SM) (B), and SM/Chol (C) lipid environments. The protein is shown in cartoon representation. The lipids are shown as sticks with N and P atoms of the head groups rendered as spheres. In (C), cholesterol aromatic rings are shown as colored planes. Orientation of the protein is different in different panels. The open form of the protein is hardly distinguishable visually and thus not shown.</figcaption></figure>
Van der Waals interactions was used in the region between 1.0 and 1.2 nm. Long-range electrostatics was computed with the Particle-Mesh-Ewald (PME) method [40] with the cut-off of explicit short-range electrostatic interactions at 1.2 nm. Velocity rescale thermostat [41] was used at the temperature of 310K. The CHARMM36 force field [42] was used for all components of the system. An integration step of 2 fs was used in all simulations with the bonds to hydrogen atoms converted to rigid constraints. Analysis was performed by the custom scripts based on the Pteros molecular modeling library [43].

## 2.5. Ensemble Docking Simulations

In order to account for protein conformational flexibility, multiple conformations from MD trajectories were extracted and used in docking simulations. The last 50 ns of MD trajectories were used to extract 834 equally spaced frames. For each frame, the docking of all ligands was performed. The docking volume was centered in the center of mass of the binding site and had dimensions 2 × 2 × 2 nm. Quick Vina 2 [44] docking software was used with default scoring function. The best docking pose was kept in each simulation. The methodology and scripts for ensemble docking used in our previous works [29–31] were used. No side chains were treated as flexible, because the possible protein flexibility was already accounted for by means of MD simulations.

# 3. Results and Discussion

## 3.1. Protein Structure and Flexibility

Both open and closed forms of the protein were found to be stable in all membrane environments (Figure 1). However, some differences are worth to be noted, as evident from the analysis of RMSD plots (Supplementary Figure S1).

<figure><img src="image_4.png" /><figcaption>**Figure 1.** Snapshots of equilibrated closed form of protein in phosphatidylcholine (PC) (A), sphingomyelin (SM) (B), and SM/Chol (C) lipid environments. The protein is shown in cartoon representation. The lipids are shown as sticks with N and P atoms of the head groups rendered as spheres. In (C), cholesterol aromatic rings are shown as colored planes. Orientation of the protein is different in different panels. The open form of the protein is hardly distinguishable visually and thus not shown.</figcaption></figure>
First, the equilibration of the protein was significantly longer for the open form than for the closed one. It took up to 400 ns for the open structure in a raft environment. In contrast, the closed form equilibrated in only 200 ns in all three membrane environments. This suggests that the closed crystal structure is closer to equilibrium conformation in a native membrane environment than the open one. Second, the final RMSD for the transmembrane region is less than ≈0.3 nm in the closed form and less than ≈0.5 nm in the open form. These are rather small values for a protein of this size. This means that the structure of the transmembrane part of the protein is in general quite rigid and is not prone to dramatic conformational changes in the course of simulations. The same trend is observed for the active site of the protein, which suggests that the dynamic of the active site is not significantly different from the rest of the protein. The protein termini (residues 1–19 and 301–342), which are exposed to aqueous solution, fluctuate randomly and interact with the membrane transiently. The interaction energies of N and C-terminal regions of the protein with the membranes are shown in Figure S2. The C-terminal domain of P2Y12 binds much stronger to the POPC membrane in both closed and open protein forms in comparison to other membrane compositions. Due to limited sampling, it is not possible to deduce whether this observation is a random coincidence or a systematic effect, but it does not correlate with the overall protein flexibility, behavior of its active site, or the binding propensity of the ligands. Thus, our simulations provide no evidence of a possible direct or indirect influence of the protein termini on the ligand binding and the flexibility of the transmembrane region of P2Y12. Despite overall stability, the protein structure is by no means “frozen” and exhibits enough local conformational mobility as it is evident from the root mean square fluctuations per residue (RMSF), which was computed for the last 50 ns of production trajectories (Figure 2).

<figure><img src="image_5.png" /><figcaption>**Figure 2.** Root means square fluctuations (RMSFs) of the transmembrane region of P2Y12 in different lipid environments for closed (A) and open (C) forms of the receptor. Differences of RMSFs between SM or raft environments with the reference PC environment for closed (B) and open (D) forms of the receptor. Filled areas correspond to the regions that are the most sensitive to the change of lipid environment.</figcaption></figure>

The pattern of root mean square fluctuations (RMSFs) in the closed form of the receptor is very similar in all lipid environments (Figure 2A). However, some notable differences
First, the equilibration of the protein was significantly longer for the open form than for the closed one. It took up to 400 ns for the open structure in a raft environment. In contrast, the closed form equilibrated in only 200 ns in all three membrane environments. This suggests that the closed crystal structure is closer to equilibrium conformation in a native membrane environment than the open one. Second, the final RMSD for the transmembrane region is less than ≈0.3 nm in the closed form and less than ≈0.5 nm in the open form. These are rather small values for a protein of this size. This means that the structure of the transmembrane part of the protein is in general quite rigid and is not prone to dramatic conformational changes in the course of simulations. The same trend is observed for the active site of the protein, which suggests that the dynamic of the active site is not significantly different from the rest of the protein. The protein termini (residues 1–19 and 301–342), which are exposed to aqueous solution, fluctuate randomly and interact with the membrane transiently. The interaction energies of N and C-terminal regions of the protein with the membranes are shown in Figure S2. The C-terminal domain of P2Y12 binds much stronger to the POPC membrane in both closed and open protein forms in comparison to other membrane compositions. Due to limited sampling, it is not possible to deduce whether this observation is a random coincidence or a systematic effect, but it does not correlate with the overall protein flexibility, behavior of its active site, or the binding propensity of the ligands. Thus, our simulations provide no evidence of a possible direct or indirect influence of the protein termini on the ligand binding and the flexibility of the transmembrane region of P2Y12. Despite overall stability, the protein structure is by no means “frozen” and exhibits enough local conformational mobility as it is evident from the root mean square fluctuations per residue (RMSF), which was computed for the last 50 ns of production trajectories (Figure 2).

<figure><img src="image_5.png" /><figcaption>**Figure 2.** Root means square fluctuations (RMSFs) of the transmembrane region of P2Y12 in different lipid environments for closed (A) and open (C) forms of the receptor. Differences of RMSFs between SM or raft environments with the reference PC environment for closed (B) and open (D) forms of the receptor. Filled areas correspond to the regions that are the most sensitive to the change of lipid environment.</figcaption></figure>

The pattern of root mean square fluctuations (RMSFs) in the closed form of the receptor is very similar in all lipid environments (Figure 2A). However, some notable differences
could be better visualized on a difference plot, where the data in the PC membrane are used as a reference. The differences of RMSFs in SM and raft membranes from the PC reference are plotted in Figure 2B. It is clearly seen that the RMSFs in both SM and raft membranes are slightly smaller than those in the PC membrane (the values are mostly negative). This trend is more pronounced for the raft mixture. This shows that the closed form is more rigid in the SM and raft environments.

The fluctuations of the open form are more interesting. This structure appears as rigid as the closed one except for two distinct regions (residues 252–269 and 220–235), which are much more flexible than the rest of the protein (Figure 2C). Interestingly, the first region 252–269 corresponds exactly to the transmembrane helix 6 (TM6). This finding corroborates well with the fact that the most noticeable difference between open and closed crystal structures is observed for TM6 [25,26] The other region 220–235 corresponds to the protruding loop on the intracellular side of the protein (Figure 3). Another interesting feature is that flexibility of these two regions in the open form of the receptor varies according to the lipid environment (Figure 2C,D). The protein is found to be the most flexible in PC membrane and significantly less flexible in SM and raft membranes. The smallest flexibility is observed in the SM membrane but not in the cholesterol enriched membrane. This correlates with much smaller RMSD and much faster equilibration of the open form of protein in the SM membrane (Figure 2A,B). The drastic decrease of flexibility in SM is more than three times more pronounced than in the case of the closed form of receptor in the regions 220–235 and 252–269 (Figure 2D). Figure 3 visualizes these regions in the aligned crystal structures of the closed and open forms of the receptor.

<figure><img src="image_6.png" /><figcaption>**Figure 3.** Aligned structures of the closed (PDB code 4PXZ:A) and open (PDB code 4NTJ:A) forms of P2Y12 receptor. Only the transmembrane core (residues 20–300) is shown for clarity. The regions that exhibit the maximal decrease of flexibility in PSM and raft membrane environments are shown in yellow for open form and in red for closed form. The region on top is residues 252–269, and the region on bottom is residues 220–235.</figcaption></figure>

We performed an in-depth analysis of the conformational mobility of TM6 by computing the RMSD of the region 252–269 with both closed and open crystal structures of P2Y12 after alignment on the rest of the membrane core of the protein (residues 20–300). Only equilibrated parts of the trajectories were used for analysis. The results are shown in Figure 4.
could be better visualized on a difference plot, where the data in the PC membrane are used as a reference. The differences of RMSFs in SM and raft membranes from the PC reference are plotted in Figure 2B. It is clearly seen that the RMSFs in both SM and raft membranes are slightly smaller than those in the PC membrane (the values are mostly negative). This trend is more pronounced for the raft mixture. This shows that the closed form is more rigid in the SM and raft environments.

The fluctuations of the open form are more interesting. This structure appears as rigid as the closed one except for two distinct regions (residues 252–269 and 220–235), which are much more flexible than the rest of the protein (Figure 2C). Interestingly, the first region 252–269 corresponds exactly to the transmembrane helix 6 (TM6). This finding corroborates well with the fact that the most noticeable difference between open and closed crystal structures is observed for TM6 [25,26] The other region 220–235 corresponds to the protruding loop on the intracellular side of the protein (Figure 3). Another interesting feature is that flexibility of these two regions in the open form of the receptor varies according to the lipid environment (Figure 2C,D). The protein is found to be the most flexible in PC membrane and significantly less flexible in SM and raft membranes. The smallest flexibility is observed in the SM membrane but not in the cholesterol enriched membrane. This correlates with much smaller RMSD and much faster equilibration of the open form of protein in the SM membrane (Figure 2A,B). The drastic decrease of flexibility in SM is more than three times more pronounced than in the case of the closed form of receptor in the regions 220–235 and 252–269 (Figure 2D). Figure 3 visualizes these regions in the aligned crystal structures of the closed and open forms of the receptor.

<figure><img src="image_6.png" /><figcaption>**Figure 3.** Aligned structures of the closed (PDB code 4PXZ:A) and open (PDB code 4NTJ:A) forms of P2Y12 receptor. Only the transmembrane core (residues 20–300) is shown for clarity. The regions that exhibit the maximal decrease of flexibility in PSM and raft membrane environments are shown in yellow for open form and in red for closed form. The region on top is residues 252–269, and the region on bottom is residues 220–235.</figcaption></figure>

We performed an in-depth analysis of the conformational mobility of TM6 by computing the RMSD of the region 252–269 with both closed and open crystal structures of P2Y12 after alignment on the rest of the membrane core of the protein (residues 20–300). Only equilibrated parts of the trajectories were used for analysis. The results are shown in Figure 4.
<figure><img src="image_7.png" /><figcaption>**Figure 4.** Root mean square deviations (RMSDs) of the TM6 region with closed and open crystal structures of P2Y12. Each point corresponds to a single simulation frame and the clouds of points are colored according to system being simulated. The regions above (white) and below (shaded) the diagonal correspond to the structures with conformations of TM6 resembling one of the closed and open crystal structures, respectively.</figcaption></figure>

The conformation of TM6 helix was stable in all simulations experiencing only slight fluctuations, which is obvious from well-localized clouds of points corresponding to each trajectory. In the simulations that started from the closed structure, the orientation of the TM6 helix evolved significantly in the cases of PC and raft environments and stabilized in very similar intermediate conformations between closed and open crystal structure (overlapping gray and blue clouds in Figure 4 are close the diagonal). In the case of the SM environment, the conformation of TM6 remained closer to initial closed structure (red cloud in Figure 4) and is distinct from the other two environments.

The behavior of the TM6 in the simulations, which started from the open structure, differs significantly. Each environment possesses its own unique conformation of the TM6, which are quite different from each other. In the case of the PC environment, the conformation of the MT6 drifted quite far away from both closed and open structures (green cloud in Figure 4). In the cases of raft and SM environments, the deviations from the crystal structures are smaller. It is interesting that both clouds corresponding to SM environments (red and pink) are almost perfectly symmetric in respect to the diagonal, and the TM6 conformation in this environment remains close to initial ones (either closed or open). This suggests that the protein is the most rigid and less prone to conformational changes in the SM environment, which is in agreement with RMSF data in Figure 2.

We did not see any spontaneous transitions between open-like and closed-like con-
formations of the TM6 in our simulations. However, we observed the multitude of its
conformations, which strongly depend on initial conditions and the lipid environment.
This allows us to speculate that the TM6 may behave as a flexible "gate", especially in a
less ordered lipid environment, which allows a large conformational mobility of P2Y12.
The time scale of these gating motions is likely to be far beyond the length of trajectories in
this study and requires extensive dedicated simulations to be detected directly.

We also analyzed the dynamics of the P2Y12 binding pocket directly by means of the
principal components analysis. The three top principal components corresponding to the
largest eigenvalues were computed for each system. All atoms of the residues contributing
to the binding pocket were used after aligning trajectory frames to the average conformation
<figure><img src="image_7.png" /><figcaption>**Figure 4.** Root mean square deviations (RMSDs) of the TM6 region with closed and open crystal structures of P2Y12. Each point corresponds to a single simulation frame and the clouds of points are colored according to system being simulated. The regions above (white) and below (shaded) the diagonal correspond to the structures with conformations of TM6 resembling one of the closed and open crystal structures, respectively.</figcaption></figure>

The conformation of TM6 helix was stable in all simulations experiencing only slight fluctuations, which is obvious from well-localized clouds of points corresponding to each trajectory. In the simulations that started from the closed structure, the orientation of the TM6 helix evolved significantly in the cases of PC and raft environments and stabilized in very similar intermediate conformations between closed and open crystal structure (overlapping gray and blue clouds in Figure 4 are close the diagonal). In the case of the SM environment, the conformation of TM6 remained closer to initial closed structure (red cloud in Figure 4) and is distinct from the other two environments.

The behavior of the TM6 in the simulations, which started from the open structure, differs significantly. Each environment possesses its own unique conformation of the TM6, which are quite different from each other. In the case of the PC environment, the conformation of the MT6 drifted quite far away from both closed and open structures (green cloud in Figure 4). In the cases of raft and SM environments, the deviations from the crystal structures are smaller. It is interesting that both clouds corresponding to SM environments (red and pink) are almost perfectly symmetric in respect to the diagonal, and the TM6 conformation in this environment remains close to initial ones (either closed or open). This suggests that the protein is the most rigid and less prone to conformational changes in the SM environment, which is in agreement with RMSF data in Figure 2.

We did not see any spontaneous transitions between open-like and closed-like con-
formations of the TM6 in our simulations. However, we observed the multitude of its
conformations, which strongly depend on initial conditions and the lipid environment.
This allows us to speculate that the TM6 may behave as a flexible "gate", especially in a
less ordered lipid environment, which allows a large conformational mobility of P2Y12.
The time scale of these gating motions is likely to be far beyond the length of trajectories in
this study and requires extensive dedicated simulations to be detected directly.

We also analyzed the dynamics of the P2Y12 binding pocket directly by means of the
principal components analysis. The three top principal components corresponding to the
largest eigenvalues were computed for each system. All atoms of the residues contributing
to the binding pocket were used after aligning trajectory frames to the average conformation
of the P2Y12 transmembrane core. The normalized overlaps between the covariance
matrices were computed between different simulated systems (Supplementary Table S2).
There is no significant overlap of the top principal components in any pair of simulated
systems—all normalized overlaps are in the range 0.2–0.33 (one means that the principal
components are identical, while zero means that they are completely different). There
is no clear correlation between the principal component overlaps and the fluidity of the
membrane environment or initial protein form. Thus, our relatively short simulations do
not allow mapping the influence of the membrane environment directly onto the motions of
the protein binding pocket. Much longer simulations, which sample characteristic motions
of the pocket better, may succeed in this.

It is possible to conclude that the decrease of membrane fluidity and increase of the
membrane ordering leads to a general decrease of the protein flexibility. This effect is
much more pronounced for the open form of the protein. Particularly, the flexibility of
TM6 (residues 252–269), which is responsible for the key structural differences between
open and closed states, decreases up to 40% in an SM environment and up to 25% in a raft
environment. The flexibility of the intracellular loop region 220–235 decreases up to 70% in
the SM environment and up to 60% in the raft environment. This points out the importance
of the lipid environment on the flexibility of this protein. It also proves that the models of
P2Y12, which do not take into account of their lipid environment, could lead to artifacts.

3.2. *Docking Results*

The ensemble docking utilized in this study allowed us to sample multiple conforma-
tions of the protein for each ligand. It takes into account changes of the protein flexibility
in each of the studied lipid environments. The distributions of obtained docking scores are
shown in Figure 5.

<figure><img src="image_8.png" /><figcaption>**Figure 5.** Distributions of the docking scores for studied ligands in different membrane environments.</figcaption></figure>

The first remarkable feature is the binding scores of TIC and its metabolites. It is
clearly seen that the active metabolite M8 shows stronger binding than pristine TIC (the
distribution is shifted to the left) in all membranes and in both protein conformations. The
other metabolite M5 shows dramatically worse binding than TIC (the distribution is shifted
to the right substantially). The difference between the peaks of TIC and M8 distributions
is ≈0.3–0.5 kcal/mol, while the difference between TIC and M5 is ≈2–2.5 kcal/mol in all
membranes. This feature corroborates well with the experimental data, which show that
of the P2Y12 transmembrane core. The normalized overlaps between the covariance
matrices were computed between different simulated systems (Supplementary Table S2).
There is no significant overlap of the top principal components in any pair of simulated
systems—all normalized overlaps are in the range 0.2–0.33 (one means that the principal
components are identical, while zero means that they are completely different). There
is no clear correlation between the principal component overlaps and the fluidity of the
membrane environment or initial protein form. Thus, our relatively short simulations do
not allow mapping the influence of the membrane environment directly onto the motions of
the protein binding pocket. Much longer simulations, which sample characteristic motions
of the pocket better, may succeed in this.

It is possible to conclude that the decrease of membrane fluidity and increase of the
membrane ordering leads to a general decrease of the protein flexibility. This effect is
much more pronounced for the open form of the protein. Particularly, the flexibility of
TM6 (residues 252–269), which is responsible for the key structural differences between
open and closed states, decreases up to 40% in an SM environment and up to 25% in a raft
environment. The flexibility of the intracellular loop region 220–235 decreases up to 70% in
the SM environment and up to 60% in the raft environment. This points out the importance
of the lipid environment on the flexibility of this protein. It also proves that the models of
P2Y12, which do not take into account of their lipid environment, could lead to artifacts.

3.2. *Docking Results*

The ensemble docking utilized in this study allowed us to sample multiple conforma-
tions of the protein for each ligand. It takes into account changes of the protein flexibility
in each of the studied lipid environments. The distributions of obtained docking scores are
shown in Figure 5.

<figure><img src="image_8.png" /><figcaption>**Figure 5.** Distributions of the docking scores for studied ligands in different membrane environments.</figcaption></figure>

The first remarkable feature is the binding scores of TIC and its metabolites. It is
clearly seen that the active metabolite M8 shows stronger binding than pristine TIC (the
distribution is shifted to the left) in all membranes and in both protein conformations. The
other metabolite M5 shows dramatically worse binding than TIC (the distribution is shifted
to the right substantially). The difference between the peaks of TIC and M8 distributions
is ≈0.3–0.5 kcal/mol, while the difference between TIC and M5 is ≈2–2.5 kcal/mol in all
membranes. This feature corroborates well with the experimental data, which show that
the primary metabolite of ticagrelor (M8) is also active with potency equal to the parent compound [12].

The specific agonist 2MeS-ADP binds a bit worse than pristine ADP in most cases; the difference between them varies from zero to ≈0.5 kcal/mol. This difference is systematic and thus unlikely to be an artifact of our docking methodology, but its relevance for the real receptor is not clear. The binding energy of the AZD1283 antagonist varies from one system to the other but remains comparable to or slightly better than that of TIC and much better than that of ADP. One notable exception is the raft membrane in closed form, where AZD1283 binds ≈0.5 kcal/mol weaker than TIC.

In general, the binding scores of TIC, its active metabolite M8, and AZD1283 are better than those of ADP, which means that these compounds could serve as concurrent antagonists. This correlates well with the experimentally observed action of these compounds. There is a strong and systematic influence of the lipid environment on the ligand binding—the binding of all ligands becomes stronger as the membrane rigidity increases. Figure 6 shows the best and average binding scores for all studied systems. The binding becomes stronger in a raw PC → SM → raft, i.e., with the increase of membrane rigidity.

<figure><img src="image_9.png" /><figcaption>**Figure 6.** Best and average binding scores for all studied systems. (A) The best binding scores. (B) The average binding scores with the corresponding standard deviations. The ligands are shown on the bottom and three different membrane environments are shown on top of each ligand. The data points are colored according to membrane environments for clarity (black for PC, red for SM and blue for the raft).</figcaption></figure>

Another systematic trend concerns the comparison of binding with respect to the two forms of the receptor. In general, the binding of all ligands to the closed form of protein (solid symbols in Figure 6) is stronger in comparison to the open one (open symbols). This
the primary metabolite of ticagrelor (M8) is also active with potency equal to the parent compound [12].

The specific agonist 2MeS-ADP binds a bit worse than pristine ADP in most cases; the difference between them varies from zero to ≈0.5 kcal/mol. This difference is systematic and thus unlikely to be an artifact of our docking methodology, but its relevance for the real receptor is not clear. The binding energy of the AZD1283 antagonist varies from one system to the other but remains comparable to or slightly better than that of TIC and much better than that of ADP. One notable exception is the raft membrane in closed form, where AZD1283 binds ≈0.5 kcal/mol weaker than TIC.

In general, the binding scores of TIC, its active metabolite M8, and AZD1283 are better than those of ADP, which means that these compounds could serve as concurrent antagonists. This correlates well with the experimentally observed action of these compounds. There is a strong and systematic influence of the lipid environment on the ligand binding—the binding of all ligands becomes stronger as the membrane rigidity increases. Figure 6 shows the best and average binding scores for all studied systems. The binding becomes stronger in a raw PC → SM → raft, i.e., with the increase of membrane rigidity.

<figure><img src="image_9.png" /><figcaption>**Figure 6.** Best and average binding scores for all studied systems. (A) The best binding scores. (B) The average binding scores with the corresponding standard deviations. The ligands are shown on the bottom and three different membrane environments are shown on top of each ligand. The data points are colored according to membrane environments for clarity (black for PC, red for SM and blue for the raft).</figcaption></figure>

Another systematic trend concerns the comparison of binding with respect to the two forms of the receptor. In general, the binding of all ligands to the closed form of protein (solid symbols in Figure 6) is stronger in comparison to the open one (open symbols). This
effect is especially pronounced for ticagrelor and its metabolites M5 and M8 and the least pronounced for AZD1283. This is in striking contrast to the results of previous docking studies, where the static crystal structures of P2Y12 were used [27]. In these studies, the docking of ticagrelor and similar compounds with bulky N6 substituents to the closed form of protein was not successful. It was suggested that the binding pocket in the closed crystal structure is too small to fit ticagrelor and its homologues because the bulky N6 substituent clashes sterically with inwardly tilted TM6.

This demonstrates that accounting for internal protein dynamics, modulated by its lipid environment, is crucial for correct assessment of the ligand binding to P2Y12.

### 3.3. Docking to The Crystal Structures vs Ensemble Docking

The striking discrepancy between our ensemble docking results and the previous rigid docking studies stimulated us to perform an in-depth comparison between these two techniques. We performed docking to the crystal structures of closed (4PXZ:A) and open (4NTJ:A) forms of P2Y12 using the same setup as for the ensemble docking. The results are shown in Table 2.

<table><caption>**Table 2.** Best docking scores for the docking to crystal structures of closed and open forms of the protein. The corresponding best docking scores from ensemble docking in the palmitoyl-oleoyl-phosphocholine (POPC) environment are shown for comparison.</caption><thead><tr><th rowspan="2">Ligand</th><th colspan="2">Closed Form</th><th colspan="2">Open Form</th></tr><tr><th>Crystal</th><th>Ensemble</th><th>Crystal</th><th>Ensemble</th></tr></thead><tbody><tr><td>TIC</td><td>-6.3</td><td>-10.5</td><td>-9.2</td><td>-9.4</td></tr><tr><td>M5</td><td>-6.8</td><td>-8.4</td><td>-7.3</td><td>-8.0</td></tr><tr><td>M8</td><td>-4.7</td><td>-10.7</td><td>-9.3</td><td>-9.9</td></tr><tr><td>2MeS-ADP</td><td>-4.9</td><td>-8.8</td><td>-7.3</td><td>-8.1</td></tr><tr><td>ADP</td><td>-7.7</td><td>-8.8</td><td>-7.6</td><td>-8.3</td></tr><tr><td>AZD12831</td><td>-4.5</td><td>-11.1</td><td>-9.3</td><td>-10.1</td></tr></tbody></table>

The best docking scores for the crystal structures differ dramatically from the results obtained in ensemble docking. In the case of crystal structures, the binding to open form is stronger for all ligands, while for ensemble docking, the trend is strictly opposite, and the binding to closed form is stronger. The pattern of binding scores for different ligands is strikingly different as well for the closed form. Particularly M8 and AZD12831 bind substantially weaker than TIC in the crystal structures. In ensemble docking, the situation is again the opposite: M8 and AZD12831 bind stronger than TIC.

In the case of the open form of protein, the differences between ensemble docking and the docking to crystal structure are less pronounced, and the trends for different ligands are much more consistent.

In the closed form of the protein, the volume of the binding pocket is restricted by the TM6 helix, which is kinked inwards and may impose steric restraints on the ligands. The static crystal structure is especially unfortunate in terms of such restraints. The ligands do not fit well into the pocket and exhibit weak binding. MD simulations allow equilibrating positions of TM6 side chains and to sample different orientations of this helix (see Figure 4), including those that are free from the steric clashes with the ligands. As a result, the ensemble docking reveals favorable docking poses with much better binding scores and demonstrates completely different trends among the ligands.

These observations emphasize the fact that collecting statistics over many protein conformations in ensemble docking allows comparing docking score distributions as a whole, which leads to reliable and unambiguous results. The dynamics of the protein, particularly its flexibility and fluctuations, are affected dramatically by the lipid environment. Such dynamical effect could never be correctly rendered by the docking to single protein conformation. Docking to the single crystal structures of P2Y12 should be avoided because it produces incorrect and misleading results, as it is evident from Table 2.
effect is especially pronounced for ticagrelor and its metabolites M5 and M8 and the least pronounced for AZD1283. This is in striking contrast to the results of previous docking studies, where the static crystal structures of P2Y12 were used [27]. In these studies, the docking of ticagrelor and similar compounds with bulky N6 substituents to the closed form of protein was not successful. It was suggested that the binding pocket in the closed crystal structure is too small to fit ticagrelor and its homologues because the bulky N6 substituent clashes sterically with inwardly tilted TM6.

This demonstrates that accounting for internal protein dynamics, modulated by its lipid environment, is crucial for correct assessment of the ligand binding to P2Y12.

### 3.3. Docking to The Crystal Structures vs Ensemble Docking

The striking discrepancy between our ensemble docking results and the previous rigid docking studies stimulated us to perform an in-depth comparison between these two techniques. We performed docking to the crystal structures of closed (4PXZ:A) and open (4NTJ:A) forms of P2Y12 using the same setup as for the ensemble docking. The results are shown in Table 2.

<table><caption>**Table 2.** Best docking scores for the docking to crystal structures of closed and open forms of the protein. The corresponding best docking scores from ensemble docking in the palmitoyl-oleoyl-phosphocholine (POPC) environment are shown for comparison.</caption><thead><tr><th rowspan="2">Ligand</th><th colspan="2">Closed Form</th><th colspan="2">Open Form</th></tr><tr><th>Crystal</th><th>Ensemble</th><th>Crystal</th><th>Ensemble</th></tr></thead><tbody><tr><td>TIC</td><td>-6.3</td><td>-10.5</td><td>-9.2</td><td>-9.4</td></tr><tr><td>M5</td><td>-6.8</td><td>-8.4</td><td>-7.3</td><td>-8.0</td></tr><tr><td>M8</td><td>-4.7</td><td>-10.7</td><td>-9.3</td><td>-9.9</td></tr><tr><td>2MeS-ADP</td><td>-4.9</td><td>-8.8</td><td>-7.3</td><td>-8.1</td></tr><tr><td>ADP</td><td>-7.7</td><td>-8.8</td><td>-7.6</td><td>-8.3</td></tr><tr><td>AZD12831</td><td>-4.5</td><td>-11.1</td><td>-9.3</td><td>-10.1</td></tr></tbody></table>

The best docking scores for the crystal structures differ dramatically from the results obtained in ensemble docking. In the case of crystal structures, the binding to open form is stronger for all ligands, while for ensemble docking, the trend is strictly opposite, and the binding to closed form is stronger. The pattern of binding scores for different ligands is strikingly different as well for the closed form. Particularly M8 and AZD12831 bind substantially weaker than TIC in the crystal structures. In ensemble docking, the situation is again the opposite: M8 and AZD12831 bind stronger than TIC.

In the case of the open form of protein, the differences between ensemble docking and the docking to crystal structure are less pronounced, and the trends for different ligands are much more consistent.

In the closed form of the protein, the volume of the binding pocket is restricted by the TM6 helix, which is kinked inwards and may impose steric restraints on the ligands. The static crystal structure is especially unfortunate in terms of such restraints. The ligands do not fit well into the pocket and exhibit weak binding. MD simulations allow equilibrating positions of TM6 side chains and to sample different orientations of this helix (see Figure 4), including those that are free from the steric clashes with the ligands. As a result, the ensemble docking reveals favorable docking poses with much better binding scores and demonstrates completely different trends among the ligands.

These observations emphasize the fact that collecting statistics over many protein conformations in ensemble docking allows comparing docking score distributions as a whole, which leads to reliable and unambiguous results. The dynamics of the protein, particularly its flexibility and fluctuations, are affected dramatically by the lipid environment. Such dynamical effect could never be correctly rendered by the docking to single protein conformation. Docking to the single crystal structures of P2Y12 should be avoided because it produces incorrect and misleading results, as it is evident from Table 2.
### 3.4. *Comparison of the Binding Sites*

Due to the large number of docking simulations, we can estimate the probability of each ligand to be in contact with particular protein residues. The results of this analysis are shown in Supplementary Figure S3.

The pattern of binding probabilities changes with lipid environment and the protein form, but these changes are surprisingly similar for different ligands despite their different chemical nature. Detailed analysis of observed fine differences between the ligands is beyond the scope of this work.

The number of residues with high binding propensity is significantly larger in a PC environment in comparison to SM and raft environments. Particularly, in the closed form of protein in the PC environment, there are two additional groups of residues with high binding probability: 111–115 and 197,201,240,248 respectively. In the open form of protein in the PC environment, there is another additional group: 154,155,157. The raft environment exhibits two residues with high binding propensity, which are absent in other environments: 88 for closed and 18 for open form, respectively.

However, these results should be interpreted with caution because they could be biased by the limited sampling in the particular MD trajectory. Ideally, multiple independent MD and ensemble docking simulations should be performed in each lipid environment, but this requires a prohibitively large amount of time and computational resources.

Comparison with crystallographic structures show that there is a partial match between the residues with highest binding probability in our simulations and the binding pockets in closed and open crystal structures, which correspond to the complexes of P2Y12 with the agonist 2MeS-ADP and the antagonist AZD1283, respectively. In the closed structure, the residues reported to bind with 2MeS-ADP are 19, 93, 97*, 101*, 102, 105,106, 109, 152, 155*, 156, 159, 163, 175*, 179, 187*, 190*, 191, 194*, 256, 259*, 263, and 280* (the matching residues in our docking simulations are marked by asterisks). In the open structure, the residues reported to bind with AZD1283 are 101*, 102, 105, 106, 109, 152, 155*, 156, 159, 187*, 190*, 191, 194*, 195, 252, 255*, 256, 259*, 260, 276*, 279*, 280*, and 283.

Spatial locations of the binding pockets for 2MeS-ADP and AZD1283 in our simulations are similar to ones present in the crystal structures, but there is no exact correspondence, as it is evident from Figure 7. Different sets of residues that are involved in ligand binding are observed in each lipid environment.

<figure><img src="image_10.png" /><figcaption>**Figure 7.** Binding pockets for 2MeS-ADP and AZD1283 in the corresponding crystal. Green corresponds to the smallest probability; red corresponds to the largest. In the case of the crystal structure, the uniform gray color is used.</figcaption></figure>
### 3.4. *Comparison of the Binding Sites*

Due to the large number of docking simulations, we can estimate the probability of each ligand to be in contact with particular protein residues. The results of this analysis are shown in Supplementary Figure S3.

The pattern of binding probabilities changes with lipid environment and the protein form, but these changes are surprisingly similar for different ligands despite their different chemical nature. Detailed analysis of observed fine differences between the ligands is beyond the scope of this work.

The number of residues with high binding propensity is significantly larger in a PC environment in comparison to SM and raft environments. Particularly, in the closed form of protein in the PC environment, there are two additional groups of residues with high binding probability: 111–115 and 197,201,240,248 respectively. In the open form of protein in the PC environment, there is another additional group: 154,155,157. The raft environment exhibits two residues with high binding propensity, which are absent in other environments: 88 for closed and 18 for open form, respectively.

However, these results should be interpreted with caution because they could be biased by the limited sampling in the particular MD trajectory. Ideally, multiple independent MD and ensemble docking simulations should be performed in each lipid environment, but this requires a prohibitively large amount of time and computational resources.

Comparison with crystallographic structures show that there is a partial match between the residues with highest binding probability in our simulations and the binding pockets in closed and open crystal structures, which correspond to the complexes of P2Y12 with the agonist 2MeS-ADP and the antagonist AZD1283, respectively. In the closed structure, the residues reported to bind with 2MeS-ADP are 19, 93, 97*, 101*, 102, 105,106, 109, 152, 155*, 156, 159, 163, 175*, 179, 187*, 190*, 191, 194*, 256, 259*, 263, and 280* (the matching residues in our docking simulations are marked by asterisks). In the open structure, the residues reported to bind with AZD1283 are 101*, 102, 105, 106, 109, 152, 155*, 156, 159, 187*, 190*, 191, 194*, 195, 252, 255*, 256, 259*, 260, 276*, 279*, 280*, and 283.

Spatial locations of the binding pockets for 2MeS-ADP and AZD1283 in our simulations are similar to ones present in the crystal structures, but there is no exact correspondence, as it is evident from Figure 7. Different sets of residues that are involved in ligand binding are observed in each lipid environment.

<figure><img src="image_10.png" /><figcaption>**Figure 7.** Binding pockets for 2MeS-ADP and AZD1283 in the corresponding crystal. Green corresponds to the smallest probability; red corresponds to the largest. In the case of the crystal structure, the uniform gray color is used.</figcaption></figure>
**4. Limitations and Perspectives**

The major limitation of the current work is limited sampling caused by rather short MD trajectories. These trajectories are sufficient to generate useful ensemble of conformations for docking, which is confirmed by the fact that we managed to dock all ligands successfully to both forms of protein. However, much longer trajectories are needed for comprehensive study of the protein dynamics in the membrane environment. A typical timescale of the physiological conformational changes in the family of GPCR proteins is likely to be in the milliseconds range, which is not accessible in all-atom MD simulations. We presume that in the present study, we mostly see the influence of the hindrance, caused by the lipid environment, on the relaxations of two initial states, which are represented by the closed and open crystal structures. Indeed, the lateral diffusion of the lipids is an order of magnitude slower in SM and raft membranes in comparison to the PC one (Supplementary Table S1), which is likely to dump volumetric relaxations of the protein in these more rigid membranes effectively. Nevertheless, even such limited sampling of the protein dynamics reveals a pronounced effect of the membrane fluidity and ordering on the ligand binding.

Thus, our work provides proof of the principle for the influence of the membrane environment on the ligand-binding propensity of the P2Y12 receptor. Further simulations with larger time scales, multiple replications, and realistic membrane compositions are needed to determine the concrete mechanism of this phenomenon.

**5. Conclusions**

The MD and ensemble docking simulations of the full-length P2Y12 receptor in different membrane environments show that the active ticagrelor metabolite M8 (AR-C124910XX) binds to the receptor stronger than ticagrelor itself, while the metabolite M5 (AR-C133913XX) shows much weaker binding. With the increase of membrane stiffness and ordering, the binding of ticagrelor and both its metabolites becomes stronger with the best binding scores observed in the raft-like sphingomyelin-cholesterol membrane. Our results suggest that this effect is caused by the decrease of the overall protein flexibility and relaxations in the rigid membrane environment rather than due to distinct conformational changes or the rearrangement of its active site. Our data provide first direct evidence of the strong influence of the membrane environment on the ligand binding in P2Y12 receptors and suggest an optimal functioning of this protein in the lipid rafts.

**Supplementary Materials:** The following are available online at https://www.mdpi.com/1999-4923/13/4/524/s1, Figure S1: Root mean square deviations of the transmembrane part of P2Y12, Figure S2: Interaction energies of the extracellular (residues 1–19) and intracellular (residues 301–342) P2Y12 termini with the lipid bilayer for different membrane compositions and protein forms. Figure S3: Probabilities of particular residues to be in contact with docked ligands in ensemble docking simulations. Table S1: Lateral diffusion coefficients of the lipids in the studied membranes and Table S2: Normalized overlaps between the covariance matrices computed for top three principal components of the P2Y12 binding pocket.

**Author Contributions:** Conceptualization, S.D. and C.R.; methodology, F.H.; S.Y. and C.R.; software, S.Y. and C.R.; validation, S.Y.; F.H. and C.R.; formal analysis, F.H. and S.Y.; data curation, F.H. and S.Y.; writing—original draft preparation, F.H.; S.Y. and C.R.; writing—review and editing, S.D. and C.R.; funding acquisition, S.Y.; C.R. and S.D. All authors have read and agreed to the published version of the manuscript.

**Funding:** S.Y. and C.R. were supported by NATO Science for Peace and Security Programmes SPS 985291. S.Y. was supported by the European Union's Horizon 2020 research and Innovation Programme under the Marie Skłodowska-Curie grant agreement No. 796245.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** Data are available on demand.
**4. Limitations and Perspectives**

The major limitation of the current work is limited sampling caused by rather short MD trajectories. These trajectories are sufficient to generate useful ensemble of conformations for docking, which is confirmed by the fact that we managed to dock all ligands successfully to both forms of protein. However, much longer trajectories are needed for comprehensive study of the protein dynamics in the membrane environment. A typical timescale of the physiological conformational changes in the family of GPCR proteins is likely to be in the milliseconds range, which is not accessible in all-atom MD simulations. We presume that in the present study, we mostly see the influence of the hindrance, caused by the lipid environment, on the relaxations of two initial states, which are represented by the closed and open crystal structures. Indeed, the lateral diffusion of the lipids is an order of magnitude slower in SM and raft membranes in comparison to the PC one (Supplementary Table S1), which is likely to dump volumetric relaxations of the protein in these more rigid membranes effectively. Nevertheless, even such limited sampling of the protein dynamics reveals a pronounced effect of the membrane fluidity and ordering on the ligand binding.

Thus, our work provides proof of the principle for the influence of the membrane environment on the ligand-binding propensity of the P2Y12 receptor. Further simulations with larger time scales, multiple replications, and realistic membrane compositions are needed to determine the concrete mechanism of this phenomenon.

**5. Conclusions**

The MD and ensemble docking simulations of the full-length P2Y12 receptor in different membrane environments show that the active ticagrelor metabolite M8 (AR-C124910XX) binds to the receptor stronger than ticagrelor itself, while the metabolite M5 (AR-C133913XX) shows much weaker binding. With the increase of membrane stiffness and ordering, the binding of ticagrelor and both its metabolites becomes stronger with the best binding scores observed in the raft-like sphingomyelin-cholesterol membrane. Our results suggest that this effect is caused by the decrease of the overall protein flexibility and relaxations in the rigid membrane environment rather than due to distinct conformational changes or the rearrangement of its active site. Our data provide first direct evidence of the strong influence of the membrane environment on the ligand binding in P2Y12 receptors and suggest an optimal functioning of this protein in the lipid rafts.

**Supplementary Materials:** The following are available online at https://www.mdpi.com/1999-4923/13/4/524/s1, Figure S1: Root mean square deviations of the transmembrane part of P2Y12, Figure S2: Interaction energies of the extracellular (residues 1–19) and intracellular (residues 301–342) P2Y12 termini with the lipid bilayer for different membrane compositions and protein forms. Figure S3: Probabilities of particular residues to be in contact with docked ligands in ensemble docking simulations. Table S1: Lateral diffusion coefficients of the lipids in the studied membranes and Table S2: Normalized overlaps between the covariance matrices computed for top three principal components of the P2Y12 binding pocket.

**Author Contributions:** Conceptualization, S.D. and C.R.; methodology, F.H.; S.Y. and C.R.; software, S.Y. and C.R.; validation, S.Y.; F.H. and C.R.; formal analysis, F.H. and S.Y.; data curation, F.H. and S.Y.; writing—original draft preparation, F.H.; S.Y. and C.R.; writing—review and editing, S.D. and C.R.; funding acquisition, S.Y.; C.R. and S.D. All authors have read and agreed to the published version of the manuscript.

**Funding:** S.Y. and C.R. were supported by NATO Science for Peace and Security Programmes SPS 985291. S.Y. was supported by the European Union's Horizon 2020 research and Innovation Programme under the Marie Skłodowska-Curie grant agreement No. 796245.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** Data are available on demand.
**Acknowledgments:** The computations were performed using the Centre de calcul Champagne-Ardenne Romeo and the Mésocentre de calcul de Franche-Comté.

**Conflicts of Interest:** The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

## References

1. Born, G.V. Adenosine Diphosphate as a Mediator of Platelet Aggregation In Vivo: An Editorial View. *Circulation* **1985**, *72*, 741–746. [CrossRef] [PubMed]

2. Saboor, M.; Ayub, Q.; Ilyas, S. Moinuddin Platelet receptors; an instrumental of platelet physiology. *Pak. J. Med. Sci.* **2013**, *29*, 891–896. [CrossRef] [PubMed]

3. Cattaneo, M. *14-The Platelet P2 Receptors*; Academic Press: Delft, The Netherlands, 2019; pp. 259–277.

4. Kunapuli, S.P.; Daniel, J.L. P2 Receptor Subtypes in the Cardiovascular System. *Biochem. J.* **1998**, *336*, 513–523. [CrossRef] [PubMed]

5. Cusack, N.J.; Hourani, S.M.O. Platelet P2 Receptors: From Curiosity to Clinical Targets. *J. Auton. Nerv. Syst.* **2000**, *81*, 37–43. [CrossRef]

6. Haberstock-Debic, H.; Andre, P.; Mills, S.; Phillips, D.R.; Conley, P.B. A Clopidogrel-Insensitive Inducible Pool of P2Y12 Receptors Contributes to Thrombus Formation: Inhibition by Elinogrel, a Direct-Acting, Reversible P2Y12 Antagonist. *J. Pharm. Exp. 2011*, *339*, 54–61. [CrossRef]

7. Savi, P.; Herbert, J.-M. Clopidogrel and Ticlopidine: P2Y12 Adenosine Diphosphate-Receptor Antagonists for the Prevention of Atherothrombosis. *Semin. Thromb. Hemost.* **2005**, *31*, 174–183. [CrossRef]

8. Angiolillo, D.J.; Fernandez-Ortiz, A.; Bernardo, E.; Alfonso, F.; Macaya, C.; Bass, T.A.; Costa, M.A. Variability in Individual Responsiveness to Clopidogrel. *J. Am. Coll. Cardiol.* **2007**, *49*, 1505–1516. [CrossRef]

9. Algaier, I.; Jakubowski, J.A.; Asai, F.; Von Kügelgen, I. Interaction of the Active Metabolite of Prasugrel, R-138727, with Cysteine 97 and Cysteine 175 of the Human P2Y12 Receptor. *J. Thromb. Haemost.* **2008**, *6*, 1908–1914. [CrossRef]

10. Savi, P.; Zachayus, J.L.; Delesque-Touchard, N.; Labouret, C.; Herve, C.; Uzabiaga, M.F.; Pereillo, J.M.; Culouscou, J.M.; Bono, F.; Ferrara, P.; et al. The Active Metabolite of Clopidogrel Disrupts P2Y12 Receptor Oligomers and Partitions Them out of Lipid Rafts. *Proc. Natl. Acad. Sci. USA* **2006**, *103*, 11069–11074. [CrossRef]

11. VAN Giezen, J.J.J.; Nilsson, L.; Berntsson, P.; Wissing, B.-M.; Giordanetto, F.; Tomlinson, W.; Greasley, P.J. Ticagrelor Binds to Human P2Y (12) Independently from ADP but Antagonizes ADP-Induced Receptor Signaling and Platelet Aggregation. *J. Thromb. Haemost.* **2009**, *7*, 1556–1565. [CrossRef]

12. Teng, R.; Oliver, S.; Hayes, M.A.; Butler, K. Absorption, Distribution, Metabolism, and Excretion of Ticagrelor in Healthy Subjects. *Drug Metab. Dispos.* **2010**, *38*, 1514–1521. [CrossRef]

13. Angiolillo, D.J.; Guzman, L.A. Clinical Overview of Promising Nonthienopyridine Antiplatelet Agents. *Am. Heart J.* **2008**, *156*, 23–28. [CrossRef]

14. Kubisa, M.J.; Jezewski, M.P.; Gasecka, A.; Siller-Matula, J.M.; Postuła, M. Ticagrelor—Toward More Efficient Platelet Inhibition and Beyond. *Ther. Clin. Risk Manag.* **2018**, *14*, 129–140. [CrossRef]

15. Wallentin, L.; Becker, R.C.; Budaj, A.; Cannon, C.P.; Emanuelsson, H.; Held, C.; Morrow, J.; Husted, S.; James, S.; Katus, H.; et al. Ticagrelor versus Clopidogrel in Patients with Acute Coronary Syndromes. *N. Engl. J. Med.* **2009**, *361*, 1045–1057. [CrossRef]

16. Abbracchio, M.P.; Burnstock, G.; Boeynaems, J.-M.; Barnard, E.A.; Boyer, J.L.; Kennedy, C.; Knight, G.E.; Fumagalli, M.; Gachet, C.; Jacobson, K.A.; et al. International Union of Pharmacology LVIII: Update on the P2Y G Protein-Coupled Nucleotide Receptors: From Molecular Mechanisms and Pathophysiology to Therapy. *Pharmacol. Rev.* **2006**, *58*, 281–341. [CrossRef]

17. Schmidt, P.; Ritscher, L.; Dong, E.N.; Hermsdorf, T.; Coster, M.; Wittkopf, D.; Meiler, J.; Schoneberg, T. Identification of Determinants Required for Agonistic and Inverse Agonistic Ligand Properties at the ADP Receptor P2Y12. *Mol. Pharmacol.* **2012**, *112*, 082198. [CrossRef]

18. Bogdanov, M.; Mileykovskaya, E.; Dowhan, W. Lipids in the Assembly of Membrane Proteins and Organization of Protein Supercomplexes: Implications for Lipid-linked Disorders. In *Lipids in Health and Disease*; Quinn, P.J., Wang, X., Eds.; Springer: Dordrecht, The Netherlands, 2008; pp. 197–239.

19. Rabani, V.; Lagoutte-Renosi, J.; Series, J.; Valot, B.; Xuereb, J.-M.; Davani, S. Cholesterol-Rich Microdomains Contribute to PAR1 Signaling in Platelets Despite a Weak Localization of the Receptor in These Microdomains. *Int. J. Mol. Sci.* **2020**, *21*, 8065. [CrossRef]

20. Lagoutte-Renosi, J.; Allemand, F.; Ramseyer, C.; Rabani, V.; Davani, S. Influence of Antiplatelet Agents on the Lipid Composition of Platelet Plasma Membrane: A Lipidomics Approach with Ticagrelor and Its Active Metabolite. *Int. J. Mol. Sci.* **2021**, *22*, 1432. [CrossRef]

21. Quinton, T.M.; Kim, S.; Jin, J.; Kunapuli, S.P. Lipid Rafts Are Required in Gαi Signaling Downstream of the P2Y12 Receptor during ADP-Mediated Platelet Activation. *J. Thromb. Haemost.* **2005**, *3*, 1036–1041. [CrossRef]

22. Haghighi, F.; Rabani, V.; Pais-De-Barros, J.-P.; Davani, S. Reorganization of Platelet Membrane Sphingomyelins by Adenosine Diphosphate and Ticagrelor. *Chem. Phys. Lipids* **2018**, *216*, 25–29. [CrossRef]
**Acknowledgments:** The computations were performed using the Centre de calcul Champagne-Ardenne Romeo and the Mésocentre de calcul de Franche-Comté.

**Conflicts of Interest:** The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

## References

1. Born, G.V. Adenosine Diphosphate as a Mediator of Platelet Aggregation In Vivo: An Editorial View. *Circulation* **1985**, *72*, 741–746. [CrossRef] [PubMed]

2. Saboor, M.; Ayub, Q.; Ilyas, S. Moinuddin Platelet receptors; an instrumental of platelet physiology. *Pak. J. Med. Sci.* **2013**, *29*, 891–896. [CrossRef] [PubMed]

3. Cattaneo, M. *14-The Platelet P2 Receptors*; Academic Press: Delft, The Netherlands, 2019; pp. 259–277.

4. Kunapuli, S.P.; Daniel, J.L. P2 Receptor Subtypes in the Cardiovascular System. *Biochem. J.* **1998**, *336*, 513–523. [CrossRef] [PubMed]

5. Cusack, N.J.; Hourani, S.M.O. Platelet P2 Receptors: From Curiosity to Clinical Targets. *J. Auton. Nerv. Syst.* **2000**, *81*, 37–43. [CrossRef]

6. Haberstock-Debic, H.; Andre, P.; Mills, S.; Phillips, D.R.; Conley, P.B. A Clopidogrel-Insensitive Inducible Pool of P2Y12 Receptors Contributes to Thrombus Formation: Inhibition by Elinogrel, a Direct-Acting, Reversible P2Y12 Antagonist. *J. Pharm. Exp. 2011*, *339*, 54–61. [CrossRef]

7. Savi, P.; Herbert, J.-M. Clopidogrel and Ticlopidine: P2Y12 Adenosine Diphosphate-Receptor Antagonists for the Prevention of Atherothrombosis. *Semin. Thromb. Hemost.* **2005**, *31*, 174–183. [CrossRef]

8. Angiolillo, D.J.; Fernandez-Ortiz, A.; Bernardo, E.; Alfonso, F.; Macaya, C.; Bass, T.A.; Costa, M.A. Variability in Individual Responsiveness to Clopidogrel. *J. Am. Coll. Cardiol.* **2007**, *49*, 1505–1516. [CrossRef]

9. Algaier, I.; Jakubowski, J.A.; Asai, F.; Von Kügelgen, I. Interaction of the Active Metabolite of Prasugrel, R-138727, with Cysteine 97 and Cysteine 175 of the Human P2Y12 Receptor. *J. Thromb. Haemost.* **2008**, *6*, 1908–1914. [CrossRef]

10. Savi, P.; Zachayus, J.L.; Delesque-Touchard, N.; Labouret, C.; Herve, C.; Uzabiaga, M.F.; Pereillo, J.M.; Culouscou, J.M.; Bono, F.; Ferrara, P.; et al. The Active Metabolite of Clopidogrel Disrupts P2Y12 Receptor Oligomers and Partitions Them out of Lipid Rafts. *Proc. Natl. Acad. Sci. USA* **2006**, *103*, 11069–11074. [CrossRef]

11. VAN Giezen, J.J.J.; Nilsson, L.; Berntsson, P.; Wissing, B.-M.; Giordanetto, F.; Tomlinson, W.; Greasley, P.J. Ticagrelor Binds to Human P2Y (12) Independently from ADP but Antagonizes ADP-Induced Receptor Signaling and Platelet Aggregation. *J. Thromb. Haemost.* **2009**, *7*, 1556–1565. [CrossRef]

12. Teng, R.; Oliver, S.; Hayes, M.A.; Butler, K. Absorption, Distribution, Metabolism, and Excretion of Ticagrelor in Healthy Subjects. *Drug Metab. Dispos.* **2010**, *38*, 1514–1521. [CrossRef]

13. Angiolillo, D.J.; Guzman, L.A. Clinical Overview of Promising Nonthienopyridine Antiplatelet Agents. *Am. Heart J.* **2008**, *156*, 23–28. [CrossRef]

14. Kubisa, M.J.; Jezewski, M.P.; Gasecka, A.; Siller-Matula, J.M.; Postuła, M. Ticagrelor—Toward More Efficient Platelet Inhibition and Beyond. *Ther. Clin. Risk Manag.* **2018**, *14*, 129–140. [CrossRef]

15. Wallentin, L.; Becker, R.C.; Budaj, A.; Cannon, C.P.; Emanuelsson, H.; Held, C.; Morrow, J.; Husted, S.; James, S.; Katus, H.; et al. Ticagrelor versus Clopidogrel in Patients with Acute Coronary Syndromes. *N. Engl. J. Med.* **2009**, *361*, 1045–1057. [CrossRef]

16. Abbracchio, M.P.; Burnstock, G.; Boeynaems, J.-M.; Barnard, E.A.; Boyer, J.L.; Kennedy, C.; Knight, G.E.; Fumagalli, M.; Gachet, C.; Jacobson, K.A.; et al. International Union of Pharmacology LVIII: Update on the P2Y G Protein-Coupled Nucleotide Receptors: From Molecular Mechanisms and Pathophysiology to Therapy. *Pharmacol. Rev.* **2006**, *58*, 281–341. [CrossRef]

17. Schmidt, P.; Ritscher, L.; Dong, E.N.; Hermsdorf, T.; Coster, M.; Wittkopf, D.; Meiler, J.; Schoneberg, T. Identification of Determinants Required for Agonistic and Inverse Agonistic Ligand Properties at the ADP Receptor P2Y12. *Mol. Pharmacol.* **2012**, *112*, 082198. [CrossRef]

18. Bogdanov, M.; Mileykovskaya, E.; Dowhan, W. Lipids in the Assembly of Membrane Proteins and Organization of Protein Supercomplexes: Implications for Lipid-linked Disorders. In *Lipids in Health and Disease*; Quinn, P.J., Wang, X., Eds.; Springer: Dordrecht, The Netherlands, 2008; pp. 197–239.

19. Rabani, V.; Lagoutte-Renosi, J.; Series, J.; Valot, B.; Xuereb, J.-M.; Davani, S. Cholesterol-Rich Microdomains Contribute to PAR1 Signaling in Platelets Despite a Weak Localization of the Receptor in These Microdomains. *Int. J. Mol. Sci.* **2020**, *21*, 8065. [CrossRef]

20. Lagoutte-Renosi, J.; Allemand, F.; Ramseyer, C.; Rabani, V.; Davani, S. Influence of Antiplatelet Agents on the Lipid Composition of Platelet Plasma Membrane: A Lipidomics Approach with Ticagrelor and Its Active Metabolite. *Int. J. Mol. Sci.* **2021**, *22*, 1432. [CrossRef]

21. Quinton, T.M.; Kim, S.; Jin, J.; Kunapuli, S.P. Lipid Rafts Are Required in Gαi Signaling Downstream of the P2Y12 Receptor during ADP-Mediated Platelet Activation. *J. Thromb. Haemost.* **2005**, *3*, 1036–1041. [CrossRef]

22. Haghighi, F.; Rabani, V.; Pais-De-Barros, J.-P.; Davani, S. Reorganization of Platelet Membrane Sphingomyelins by Adenosine Diphosphate and Ticagrelor. *Chem. Phys. Lipids* **2018**, *216*, 25–29. [CrossRef]
23. Deflorian, F.; Jacobson, K.A. Comparison of Three GPCR Structural Templates for Modeling of the P2Y12 Nucleotide Receptor. *J. Comput. Aided Mol. Des.* **2011**, *25*, 329–338. [CrossRef] [PubMed]

24. Liu, H.; Ge, H.; Peng, Y.; Xiao, P.; Xu, J. Molecular mechanism of action for reversible P2Y12 antagonists. *Biophys. Chem.* **2011**, *155*, 74–81. [CrossRef]

25. Zhang, K.; Zhang, J.; Gao, Z.-G.; Zhang, D.; Zhu, L.; Han, G.W.; Moss, S.M.; Paoletta, S.; Kiselev, E.; Lu, W.; et al. Structure of the Human P2Y12 Receptor in Complex with an Antithrombotic Drug. *Nature* **2014**, *509*, 115–118. [CrossRef] [PubMed]

26. Zhang, J.; Zhang, K.; Gao, Z.-G.; Paoletta, S.; Zhang, D.; Han, G.W.; Li, T.; Ma, L.; Zhang, W.; Müller, C.E.; et al. Agonist-Bound Structure of the Human P2Y12 Receptor. *Nature* **2014**, *509*, 119–122. [CrossRef] [PubMed]

27. Paoletta, S.; Sabbadin, D.; Kügelgen, I.; Hinz, S.; Katritch, V.; Hoffmann, K.; Abdelrahman, A.; Straßburger, J.; Baqi, Y.; Zhao, Q.; et al. Modeling Ligand Recognition at the P2Y12 Receptor in Light of X-Ray Structural Information. *J. Comput. Aided Mol. Des.* **2015**, *29*, 737–756. [CrossRef] [PubMed]

28. Yesylevskyy, S.O.; Ramseyer, C.; Pudlo, M.; Pallandre, J.-R.; Borg, C. Selective Inhibition of STAT3 with Respect to STAT1: Insights from Molecular Dynamics and Ensemble Docking Simulations. *J. Chem. Inf. Modeling* **2016**, *56*, 1588–1596. [CrossRef] [PubMed]

29. Pallandre, J.-R.; Borg, C.; Rognan, D.; Boibessot, T.; Luzet, V.; Yesylevskyy, S.; Ramseyer, C.; Pudlo, M. Novel Aminotetrazole Derivatives as Selective STAT3 Non-Peptide Inhibitors. *Eur. J. Med. Chem.* **2015**, *103*, 163–174. [CrossRef]

30. Yesylevskyy, S.O.; Ramseyer, C.; Savenko, M.; Mura, S.; Couvreur, P. Low-Density Lipoproteins and Human Serum Albumin as Carriers of Squalenoylated Drugs: Insights from Molecular Simulations. *Mol. Pharm.* **2018**, *15*, 585–591. [CrossRef]

31. Gobeaux, F.; Bizeau, J.; Samson, F.; Marichal, L.; Grillo, I.; Wien, F.; Yesylevsky, S.O.; Ramseyer, C.; Rouquette, M.; Lepêtre-Mouelhi, S.; et al. Albumin-driven disassembly of lipidic nanoparticles: The specific case of the squalene-adenosine nanodrug. *Nanoscale* **2020**, *12*, 2793–2809. [CrossRef]

32. Ko, J.; Park, H.; Heo, L.; Seok, C. GalaxyWEB Server for Protein Structure Prediction and Refinement. *Nucleic Acids Res.* **2012**, *40*, 294–297. [CrossRef]

33. Zhang, C.; Mortuza, S.M.; He, B.; Wang, Y.; Zhang, Y. Template-Based and Free Modeling of I-TASSER and QUARK Pipelines Using Predicted Contact Maps in CASP12. *Proteins Struct. Funct. Bioinform.* **2018**, *86*, 136–151. [CrossRef]

34. Lagoutte-Renosi, J.; Royer, B.; Rabani, V.; Davani, S. Validation of an HPLC-MS/MS Method for the Determination of Plasma Ticagrelor and Its Active Metabolite Useful for Research and Clinical Practice. *Molecules* **2021**, *26*, 278. [CrossRef]

35. Lee, J.; Cheng, X.; Swails, J.M.; Yeom, M.S.; Eastman, P.K.; Lemkul, J.A.; Wei, S.; Buckner, J.; Jeong, J.C.; Qi, Y.; et al. CHARMM-GUI Input Generator for NAMD, GROMACS, AMBER, OpenMM, and CHARMM/OpenMM Simulations Using the CHARMM36 Additive Force Field. *J. Chem. Theory Comput.* **2016**, *12*, 405–413. [CrossRef]

36. Jo, S.; Kim, T.; Iyer, V.G.; Im, W. CHARMM-GUI: A Web-Based Graphical User Interface for CHARMM. *J. Comput. Chem.* **2008**, *29*, 1859–1865. [CrossRef]

37. Hess, B.; Kutzner, C.; Spoel, D.; Lindahl, E. GROMACS 4: Algorithms for Highly Efficient, Load-Balanced, and Scalable Molecular Simulation. *J. Chem. Theory Comput.* **2008**, *4*, 435–447. [CrossRef]

38. Berendsen, H.J.C.; Postma, J.P.M.; Gunsteren, W.F.; DiNola, A.; Haak, J.R. Molecular Dynamics with Coupling to an External Bath. *J. Chem. Phys.* **1984**, *81*, 3684–3690. [CrossRef]

39. Páll, S.; Hess, B. A Flexible Algorithm for Calculating Pair Interactions on SIMD Architectures. *Comput. Phys. Commun.* **2013**, *184*, 2641–2650. [CrossRef]

40. Spoel, D.; Lindahl, E.; Hess, B.; Groenhof, G.; Mark, A.E.; Berendsen, H.J.C. GROMACS: Fast, Flexible and Free. *J. Comp. Chem* **2005**, *26*, 1701–1718. [CrossRef]

41. Bussi, G.; Donadio, D.; Parrinello, M. Canonical Sampling through Velocity Rescaling. *J. Chem. Phys.* **2007**, *126*, 014101. [CrossRef]

42. Huang, J.; MacKerell, A.D., Jr. CHARMM36 All-Atom Additive Protein Force Field: Validation Based on Comparison to NMR Data. *J. Comput. Chem.* **2013**, *34*, 2135–2145. [CrossRef]

43. Yesylevskyy, S.O. Pteros 2.0: Evolution of the Fast Parallel Molecular Analysis Library for C++ and Python. *J. Comput. Chem.* **2015**, *36*, 1480–1488. [CrossRef] [PubMed]

44. Alhossary, A.; Handoko, S.D.; Mu, Y.; Kwoh, C.-K. Fast, Accurate, and Reliable Molecular Docking with QuickVina 2. *Bioinformatics* **2015**, *31*, 2214–2216. [CrossRef] [PubMed]
23. Deflorian, F.; Jacobson, K.A. Comparison of Three GPCR Structural Templates for Modeling of the P2Y12 Nucleotide Receptor. *J. Comput. Aided Mol. Des.* **2011**, *25*, 329–338. [CrossRef] [PubMed]

24. Liu, H.; Ge, H.; Peng, Y.; Xiao, P.; Xu, J. Molecular mechanism of action for reversible P2Y12 antagonists. *Biophys. Chem.* **2011**, *155*, 74–81. [CrossRef]

25. Zhang, K.; Zhang, J.; Gao, Z.-G.; Zhang, D.; Zhu, L.; Han, G.W.; Moss, S.M.; Paoletta, S.; Kiselev, E.; Lu, W.; et al. Structure of the Human P2Y12 Receptor in Complex with an Antithrombotic Drug. *Nature* **2014**, *509*, 115–118. [CrossRef] [PubMed]

26. Zhang, J.; Zhang, K.; Gao, Z.-G.; Paoletta, S.; Zhang, D.; Han, G.W.; Li, T.; Ma, L.; Zhang, W.; Müller, C.E.; et al. Agonist-Bound Structure of the Human P2Y12 Receptor. *Nature* **2014**, *509*, 119–122. [CrossRef] [PubMed]

27. Paoletta, S.; Sabbadin, D.; Kügelgen, I.; Hinz, S.; Katritch, V.; Hoffmann, K.; Abdelrahman, A.; Straßburger, J.; Baqi, Y.; Zhao, Q.; et al. Modeling Ligand Recognition at the P2Y12 Receptor in Light of X-Ray Structural Information. *J. Comput. Aided Mol. Des.* **2015**, *29*, 737–756. [CrossRef] [PubMed]

28. Yesylevskyy, S.O.; Ramseyer, C.; Pudlo, M.; Pallandre, J.-R.; Borg, C. Selective Inhibition of STAT3 with Respect to STAT1: Insights from Molecular Dynamics and Ensemble Docking Simulations. *J. Chem. Inf. Modeling* **2016**, *56*, 1588–1596. [CrossRef] [PubMed]

29. Pallandre, J.-R.; Borg, C.; Rognan, D.; Boibessot, T.; Luzet, V.; Yesylevskyy, S.; Ramseyer, C.; Pudlo, M. Novel Aminotetrazole Derivatives as Selective STAT3 Non-Peptide Inhibitors. *Eur. J. Med. Chem.* **2015**, *103*, 163–174. [CrossRef]

30. Yesylevskyy, S.O.; Ramseyer, C.; Savenko, M.; Mura, S.; Couvreur, P. Low-Density Lipoproteins and Human Serum Albumin as Carriers of Squalenoylated Drugs: Insights from Molecular Simulations. *Mol. Pharm.* **2018**, *15*, 585–591. [CrossRef]

31. Gobeaux, F.; Bizeau, J.; Samson, F.; Marichal, L.; Grillo, I.; Wien, F.; Yesylevsky, S.O.; Ramseyer, C.; Rouquette, M.; Lepêtre-Mouelhi, S.; et al. Albumin-driven disassembly of lipidic nanoparticles: The specific case of the squalene-adenosine nanodrug. *Nanoscale* **2020**, *12*, 2793–2809. [CrossRef]

32. Ko, J.; Park, H.; Heo, L.; Seok, C. GalaxyWEB Server for Protein Structure Prediction and Refinement. *Nucleic Acids Res.* **2012**, *40*, 294–297. [CrossRef]

33. Zhang, C.; Mortuza, S.M.; He, B.; Wang, Y.; Zhang, Y. Template-Based and Free Modeling of I-TASSER and QUARK Pipelines Using Predicted Contact Maps in CASP12. *Proteins Struct. Funct. Bioinform.* **2018**, *86*, 136–151. [CrossRef]

34. Lagoutte-Renosi, J.; Royer, B.; Rabani, V.; Davani, S. Validation of an HPLC-MS/MS Method for the Determination of Plasma Ticagrelor and Its Active Metabolite Useful for Research and Clinical Practice. *Molecules* **2021**, *26*, 278. [CrossRef]

35. Lee, J.; Cheng, X.; Swails, J.M.; Yeom, M.S.; Eastman, P.K.; Lemkul, J.A.; Wei, S.; Buckner, J.; Jeong, J.C.; Qi, Y.; et al. CHARMM-GUI Input Generator for NAMD, GROMACS, AMBER, OpenMM, and CHARMM/OpenMM Simulations Using the CHARMM36 Additive Force Field. *J. Chem. Theory Comput.* **2016**, *12*, 405–413. [CrossRef]

36. Jo, S.; Kim, T.; Iyer, V.G.; Im, W. CHARMM-GUI: A Web-Based Graphical User Interface for CHARMM. *J. Comput. Chem.* **2008**, *29*, 1859–1865. [CrossRef]

37. Hess, B.; Kutzner, C.; Spoel, D.; Lindahl, E. GROMACS 4: Algorithms for Highly Efficient, Load-Balanced, and Scalable Molecular Simulation. *J. Chem. Theory Comput.* **2008**, *4*, 435–447. [CrossRef]

38. Berendsen, H.J.C.; Postma, J.P.M.; Gunsteren, W.F.; DiNola, A.; Haak, J.R. Molecular Dynamics with Coupling to an External Bath. *J. Chem. Phys.* **1984**, *81*, 3684–3690. [CrossRef]

39. Páll, S.; Hess, B. A Flexible Algorithm for Calculating Pair Interactions on SIMD Architectures. *Comput. Phys. Commun.* **2013**, *184*, 2641–2650. [CrossRef]

40. Spoel, D.; Lindahl, E.; Hess, B.; Groenhof, G.; Mark, A.E.; Berendsen, H.J.C. GROMACS: Fast, Flexible and Free. *J. Comp. Chem* **2005**, *26*, 1701–1718. [CrossRef]

41. Bussi, G.; Donadio, D.; Parrinello, M. Canonical Sampling through Velocity Rescaling. *J. Chem. Phys.* **2007**, *126*, 014101. [CrossRef]

42. Huang, J.; MacKerell, A.D., Jr. CHARMM36 All-Atom Additive Protein Force Field: Validation Based on Comparison to NMR Data. *J. Comput. Chem.* **2013**, *34*, 2135–2145. [CrossRef]

43. Yesylevskyy, S.O. Pteros 2.0: Evolution of the Fast Parallel Molecular Analysis Library for C++ and Python. *J. Comput. Chem.* **2015**, *36*, 1480–1488. [CrossRef] [PubMed]

44. Alhossary, A.; Handoko, S.D.; Mu, Y.; Kwoh, C.-K. Fast, Accurate, and Reliable Molecular Docking with QuickVina 2. *Bioinformatics* **2015**, *31*, 2214–2216. [CrossRef] [PubMed]