Université De Strasbourg

Iuri Casciuc

Thèse Dirigée

M Varnek

M Lehn

M Morelli

M Villoutreix

Mme Camproux

Discipline/ Spécialité : Chimie / Chémoinformatique Cartographie chimique et modélisation de systèmes complexes

Dedication

In the memory of my father Casciuc Valeriy, who passed away on the 24 th of November 2018.

Abstract

This work concerns application of Generative Topographic Mapping method to different tasks including data analysis and visualization, virtual screening and library design.

Performance of multi-target GTM-based classification models (uGTM) in virtual screening was investigated and consensus usage of several uGTMs has been suggested. Virtual screening involving a combination of GTM with some other chemoinformatics techniques allowed to discover 29 new BRD4 inhibitors, activities of which were experimentally confirmed. As a library design tool, GTM was compared to the MaxMin method. Although diversity of MaxMin libraries is systematically larger than those obtained with GTM, the latter is much faster and, therefore, can be recommended for large datasets. A modeling workflow for speciation analysis in imine-based Dynamic Combinatorial Libraries in absence and presence of a protein has been suggested. Developed models are publicly available at the site of the Laboratory of Chemoinformatics.

I would like to express my sincere gratitude to all my colleagues from the Laboratory of Chemoinformatics at the University of Strasbourg. Particular thanks to my supervisors Professor Alexandre Varnek and Professor Jean-Marie Lehn, for their patience, a plethora of professional advice and also for their kindness. Special thanks to Dr. Gilles Marcou for our talks related to the research and educational sides of the work. Special mulțumesc to Dr.

Dragos Horvath for all the tools that he has implemented, for his help with the scripts and last but not least for not only helping me to not forget the hardly-learned Romanian language but also for teaching me some new words. Dr. Igor Baskin, thank you for your advice and for our productive discussions while I was working on a diverse-library project.

I appreciate the help of Dr. Fanny Bonachera and Dr. Olga Klimchuk in organizing my working process and documents and for their help in some bureaucracy unrelated to work. I am grateful to Dr. Arkadii Lin, Dr. Alexey Orlov, Yuliana Zabolotna and William Bort for our friendship and for all the discussions that we had during this time. I thank all the former PhD students from our lab: Dr. Pavel Sidorov, Dr. Timur Gimadiev and Dr. Marta Glavatskikh. I cannot find the words to express how deeply grateful I am to all the members of the lab for their support after my father's loss. Finally, I would like to thank all members of Lehn's laboratory, especially Dr. Artem Osypenko and Bohgdan Kozybroda, for our fruitful discussions on the "purely-chemical" topics. It was worth mentioning another person: Dr. Julien Diharce, a friend of mine who encouraged and supported me from my first year at the University of Nice. Last but not least, I would like to thank the Région Grand Est for the PhD fellowship.

Résumé en français 1.Introduction

Actuellement les bases de données chimiques incluent des millions de structures de composés chimiques [1,2]. Grâce à la synthèse combinatoire et aux réacteurs en flux continu ce nombre augmente exponentiellement. Néanmoins ces chiffres sont « négligeables » en comparaison du nombre de composés que contiendrait l'espace chimique même en se limitant aux molécules d'intérêt thérapeutique, celui-ci étant estimé à 10 33 [3]. L'exploration et l'analyse de cet espace permet aux chimistes de mieux comprendre les relations structure-activité ; de plus, grâce l'analyse des régions inexplorées de l'espace chimique facilite l'innovation, en particulier pour la recherche de nouveaux candidats médicaments.

L'une des approches qui permet d'effectuer une telle analyse est la méthode de cartographie topographique générative (Generative Topographic Mapping, ou GTM) [START_REF] Bishop | GTM: The generative topographic mapping[END_REF].

Elle localise les structures chimiques, représentées par un espace de descripteur multidimensionnel initial sur un espace bidimensionnel plat, appelé une carte (Figure 1.1).

La GTM établi en fait une correspondance entre une distribution de probabilité dans l'espace initial avec une distribution bidimensionnelle sur la carte. Cette dernière est quantifiée en des points spécifiques de la carte, appelés « noeuds ». À une molécule localisée dans l'espace initial correspond une distribution de probabilité sur la carte. Les coordonnées de la molécule sur la carte correspondent au centre de gravité de la distribution de probabilité qui représente le composé sur la carte. Inversement, à chaque lieu de la carte correspond une distribution de l'espace initial qui représente une population de structures chimiques. Précédemment, cette approche a été appliquée avec succès pour visualiser et analyser des données chimiques [START_REF] Gaspar | Chemical data visualization and analysis with incremental generative topographic mapping: big data challenge[END_REF] ainsi que pour préparer des modèles prédictifs de régression ou de classification [START_REF] Gaspar | GTM-Based QSAR Models and Their Applicability Domains[END_REF].

L'objectif de cette thèse est d'explorer l'application de la méthode GTM à plusieurs tâches de la chémoinformatique : le criblage virtuel, l'analyse de l'espace chimique de systèmes complexes et la constitution d'une bibliothèque de composés chimiquement divers. La thèse est divisée en 6 Chapitres. Les chapitres 1 et 2 sont l'introduction et les méthodes, respectivement. Le chapitre 3 résume les résultats de l'utilisation simultanée de plusieurs cartes GTM dans le criblage virtuel de différentes cibles biologiques à partir de données issues de la base de données ChEMBL. Le chapitre 4 décrit le projet dédié à la conception assisté par ordinateur visant à trouver de nouveaux inhibiteurs de Bromodomaine. Le chapitre 5 est est dédié à la modélisation de la spéciation de bibliothèques combinatoires dynamiques [START_REF] Lehn | Dynamic combinatorial chemistry and virtual combinatorial libraries[END_REF]8] en absence ou en présence d'un effecteur (protéine ou métal). Ici, la préparation de modèles prédictifs sont décrits, pour le calcul du logarithme de la constante d'équilibre logKeq pour les réactions de formation d'imines, et pour la formation de complexes entre des effecteurs variés et des molécules organiques.

Dans le même chapitre la constitution d'un jeu de données structurellement divers à partir d'une base de données contenant plus de 42000 composés à l'aide d'une GTM, ainsi que la comparaison de sa performance avec la méthode traditionnelle MaxMin [START_REF] Holliday | Definitions of" dissimilarity" for dissimilaritybased compound selection[END_REF] est décrite.

Dans le chapitre 6 les conclusions générales ainsi que les perspectives sont décrites.

Résultats et discussions 1.2.1 Application du modèle consensus GTM au criblage virtuel

Il a été démontré que la méthode GTM peut être utilisée pour préparer des modèles prédictifs de régression et de classification [START_REF] Gaspar | GTM-Based QSAR Models and Their Applicability Domains[END_REF]. Dans ce projet, nous montrons que le consensus de cartes construites à partir de différents descripteurs moléculaires mais sur un même jeu de données, fournis des prédictions plus fiables par rapport à l'utilisation d'une carte unique basée sur un seul ensemble de descripteurs moléculaires. Nos études ont été menées sur les « cartes universelles » [START_REF] Sidorov | Mappability of drug-like space: towards a polypharmacologically competent map of drug-relevant compounds[END_REF] qui sont capables de distinguer les composés actifs de composés inactifs pour plus de 600 cibles biologiques, simultanément [2]. Huit cartes universelles ont été obtenues à partir de différents espaces de descripteurs utilisés.

Chaque carte universelle est capable à prédire l'activité de ligands pour plusieurs cibles avec de bonnes performances. Néanmoins des performances de prédiction pour une même cible divergent considérablement d'une carte à l'autre. Ces performances sont estimées par un paramètre statistique appelé « précision balancée » (Balanced Accuracy, ou BA) qui varie entre 0.5 (si les prédictions sont aléatoires) et 1 (si les prédictions correspondent à la l'expérience). Il a été trouvé qu'aucune des 8 cartes n'est capable de séparer à elle seule les composés actifs et inactifs pour toutes les cibles avec une BA suffisamment élevée.

Toutefois, les prédictions de chaque carte peuvent être combinées en un consensus. Ainsi, il a été observé qu'un consensus de 7 cartes est suffisant pour prédire l'activité des molécules avec une haute fiabilité (BA>0.75) sur plus de 85% des cibles, simultanément (Figure 1.

2).

Ces cartes sont complémentaires, car les cibles moins bien prédites par une carte sont mieux prédites sur une autre, en exploitant les points de vus alternatifs que représentent les descripteurs moléculaires utilisés pour chacune.

Conception assistée par ordinateur de nouveaux inhibiteurs de Bromodomaine

L'étude par criblage virtuel (Virtual Screening, ou VS) décrite ici visait à identifier de nouveaux ligands de Bromodomaine BRD4. Elle s'est appuyée sur le contenu de bases de données publiques (ChEMBL, REAXYS) pour établir, dans un premier temps, un modèle prédictif de l'activité BRD puis, dans un second temps, l'utilisation de ces modèles pour la sélection de ligands putatifs. Différentes approches chémoinformatiques (SVM, pharmacophores, GTM, docking) ont donc été utilisées pour filtrer la collection de 2 millions de composés de la société Enamine. Ce partenaire industriel a ensuite testé expérimentalement un sous-ensemble de 2992 molécules de cette sélection.

Modélisation in silico des bibliothèques combinatoires dynamiques des imines 1.2.3.1 Application de la méthode GTM à la constitution d'une bibliothèque de composés chimiquement diverse

Les bibliothèques de composés chimiquement diverses sont particulièrement importantes pour recherche en chimie médicinale car un tel ensemble de composés vise à tester le plus grand nombre d'hypothèses quand une nouvelle activité biologique est recherchée et qu'il existe peu de connaissances à priori pour faire des choix rationnels.

Habituellement, une telle collection est obtenue à partir d'une sélection des composés les plus différents les uns des autres, à partir d'une grande chimiothèque de composés accessibles commercialement ou par voie de synthèse. La dissimilarité entre deux structures chimiques est assimilée à la distance séparant celles-ci dans l'espace chimique, c'est-à-dire la distance entre les deux vecteurs de descripteurs moléculaires qui les représentent.

Trois jeux de données ont été utilisés dans ce travail : l'un contenant 154 amines, l'autre -277 aldéhydes et le troisième contenant 42658 imines qui sont les hypothétiques produits de réaction entre chacune des amines avec chacun des aldéhydes. Le but est de Comme la méthode basée sur la GTM est une méthode bien plus rapide que MaxMin, elle offre des perspectives intéressantes.

Modélisation des équilibres dans une bibliothèque combinatoire dynamique

Une bibliothèque de composés qui peuvent interagir de manière réversible les uns avec les autres s'appelle une bibliothèque combinatoire dynamique (Dynamic

Conclusions

La méthode GTM a été utilisée avec succès pour l'analyse de l'espace chimique, le criblage virtuel et la constitution d'une chimiothèque diverse. Nous avons démontré qu'un consensus de modèles individuels GTM fournis des prédictions plus fiables que chaque

Liste des presentations

Liste des publications 2 Introduction

Nowadays, chemical databases include millions of structures of chemical compounds [1,2]. Thanks to combinatorial synthesis and continuous flow reactors, this number increases exponentially. However, these numbers are "negligible" in comparison to the number of compounds that the chemical space would contain even if it were limited to molecules of therapeutic interest, this being estimated at 10 33 [3]. The exploration and analysis of chemical space allow chemists to understand structure-activity relationships better; moreover, the study of the unexplored regions of the chemical space facilitates innovation, in particular for the research of new drug candidates.

When it comes to chemical data visualization, analysis and modeling, the applied methods could be either descriptor-based or graph-based. In a descriptor-based approach, a compound is represented as a vector of descriptors, and each descriptor is describing the molecule in terms of physical or chemical properties (molecular weight, logP) and/or purely structural (number of atoms, types of bonds). These vectors of descriptors are serving as input for machine-learning algorithms. Some of these algorithms, like dimensionality reduction techniques, are designed specifically for the visualization and modeling of multidimensional data. Multiple dimensionality reduction techniques are reported in the literature: Principal Component Analysis (PCA [START_REF] Pearson | LIII. On lines and planes of closest fit to systems of points in space[END_REF][START_REF] Hotelling | Analysis of a complex of statistical variables into principal components[END_REF][START_REF] Akella | Cheminformatics approaches to analyze diversity in compound screening libraries[END_REF]), Multi-Dimensional Scaling (MDS [START_REF] Buja | Interface Foundation of America[END_REF]), Sammon mapping, Self-Organizing Maps (SOM) [START_REF] Kohonen | The self-organizing map[END_REF] due to their efficiency.

In 2001 T. Oprea [START_REF] Oprea | Chemography: The art of navigating in chemical space[END_REF] proposed a new term -"… chemography, by analogy with geography, as the art of navigating chemical space." In other words, Oprea suggested the usage of maps for navigation in chemical space. Such maps of chemical space, by analogy with the world map, should posses a universal character, i.e. the compounds are defining the "contours of the continents" while their properties will be defining the colors. The above-mentioned dimensionality reduction methods are efficient, nevertheless they are not perfect. For example, SOM is producing a 2D map which is based on a non-linear model, PCA is also able to produce a 2D map if two principal components will be taken. However, PCA is efficient with datasets having internal linear correlations [START_REF] Balakin | Pharmaceutical Data Mining[END_REF], but it could fail while representing vast multidimensional data [START_REF] Maniyar | Data Visualization during the Early Stages of Drug Discovery[END_REF]. MDS is another dimensionality reduction technique that is also linear that is using Euclidean distances [START_REF] Neal | Pattern Recognition and Machine Learning[END_REF]. Sammon maps do not allow the addition of the new data on already existing map, forcing the user to rebuild the map if the new compounds should be added [START_REF] Hutchison | Intelligent Data Engineering and Automated Learning 23[END_REF].

On the other hand, graph-based approaches represent a molecule as a graph, the atoms and bonds corresponding to the graph's nodes and edges, respectively. One way to work in graph-based chemical space is to rely on the concept of a scaffold that is defined as the "core part" of the structure with all the terminal chains removed [23]. These could be regrouped in a so-called hierarchical scaffold tree, which allows the data visualization and modeling [START_REF] Schuffenhauer | The Scaffold Tree -Visualization of the Scaffold Universe by Hierarchical Scaffold Classification[END_REF].

Generative Topographic Mapping (GTM) [START_REF] Bishop | GTM: The generative topographic mapping[END_REF] is a probabilistic extension of SOM that considers the likelihood of the training data as the objective function. Moreover, unlike to SOM, the object is not associated with one particular node. In GTM, the objects are represented as a probability distribution over all the nodes, thus creating a vector of probabilities for each data point. This vector of probabilities is used for data visualization as well as for the building classification and regression models. GTM is a versatile method that can be applied in different everyday tasks of chemoinformatics like data visualization, libraries comparison, QSAR, de novo design; therefore, GTM could be compared to a Swiss army knife.

GTM has been used in several projects for data visualization. The hierarchical GTM algorithm [START_REF] Tino | Hierarchical GTM} constructing localized nonlinear projection manifolds in a principled way[END_REF] has been used [START_REF] Schuffenhauer | The Scaffold Tree -Visualization of the Scaffold Universe by Hierarchical Scaffold Classification[END_REF] for the visualization of the active/inactive classes distribution in five different datasets issued from high-throughput screening. Large datasets of compounds (2.2M compounds) have been visualized for the first time by applying the incremental GTM [START_REF] Gaspar | Chemical data visualization and analysis with incremental generative topographic mapping: big data challenge[END_REF]. In the same work, a comparison of libraries has been made. Each library was considered as a single object of cumulated responsibilities or properties. A

Responsibility Pattern (RP) term has been introduced by Klimenko et al. [START_REF] Klimenko | Chemical space mapping and structure--activity analysis of the ChEMBL antiviral compound set[END_REF]; RP allowed to automatically detect and extract the compounds similar in the latent space. The concept of "privileged substructures" (PSM) was initially introduced by BE. Evans et al. [START_REF] Evans | Methods for drug discovery: development of potent, selective, orally effective cholecystokinin antagonists[END_REF],

referring to core structures that are recurrent in compounds active against a given target family and, therefore, associated with that biological activity. This approach has been applied in the analysis and modeling of antimalarial compounds [START_REF] Sidorov | QSAR modeling and chemical space analysis of antimalarial compounds[END_REF]. PSM has been modified by applying the retrosynthetic rules (RECAP) [START_REF] Kayastha | Privileged structural motif detection and analysis using generative topographic maps[END_REF]. The authors tried to extract the "frequent" RECAP cores to identify PSMs for inhibitors of protease, kinase and GPCRs. In recent work, GTM has been applied to visualization, analysis and comparison of the compounds tested against virus species, representatives of the Coronaviridae family [START_REF] Horvath | A Chemographic Audit of Anti-Coronavirus Structure-Activity Information from Public Databases[END_REF].

GTM has been successfully applied for QSAR and QSPR modeling. Kireeva et al.

[31] used GTM-based classification models for the prediction of the melting point of ionic liquids. Gaspar et al. have used GTM-based regression models [START_REF] Gaspar | GTM-Based QSAR Models and Their Applicability Domains[END_REF] for the modeling of stability constants for metal binders, the activity of thrombin inhibitors and aqueous solubility. In multiple works, GTM has been compared to other popular machine-learning methods: SOM [START_REF] Erwin | Self-organizing maps: ordering, convergence properties and energy functions[END_REF], Random Forest [START_REF] Liaw | Classification and regression by randomForest[END_REF], Partial Least Squares [START_REF] Wold | {PLS}-regression: a basic tool of chemometrics[END_REF], M5P regression tree [START_REF] Quinlan | Learning with continuous classes[END_REF], SVM [START_REF] Cortes | Support-vector networks[END_REF]. It has been shown across many projects that GTM is a method that can compete in terms of the performance of produces target-and property-specific models with other machine-learning methods. In their work, Sidorov et al. [START_REF] Sidorov | Mappability of drug-like space: towards a polypharmacologically competent map of drug-relevant compounds[END_REF] have shown that GTM can be successfully used as a multi-target predictive model. For instance, a dataset of 1.3M

ChEMBL compounds (version 20) corresponding to 410 targets have been covered, and approximatively 80% of these targets have been predicted with relatively high Balanced Accuracy (> 0.7). Lately, Lin et al. [START_REF] Lin | Multi-task generative topographic mapping in virtual screening[END_REF] have applied the same protocol on the ChEMBL 23 dataset and benchmarked the obtained model with popular machine-learning methods.

Conformational sampling plays an important role in medicinal chemistry. Horvath et al. have described the application of GTM to conformational sampling [START_REF] Horvath | Generative Topographic Mapping of Conformational Space[END_REF]. In this work, a set of conformers with previously calculated total energies has been used (calculations were done using the general AMBER force field [START_REF] Wang | Development and testing of a general amber force field[END_REF]). Torsion angles and some non-bonded contact energies have been used as descriptors to describe the conformers interaction fingerprints. The obtained map can be used for visualization and analysis of the "training" conformational space as well as to predict the energies for new conformers. This approach has been applied to the conformational space of dipeptides [START_REF] Horvath | Monitoring of the Conformational Space of Dipeptides by Generative Topographic Mapping[END_REF]. It was also used in a docking study of the ATP-binding site of CDK2 [START_REF] Horvath | Generative Topographic Mapping of the Docking Conformational Space[END_REF]. The maps have been trained to be able to discriminate native from non-native ligand poses as well as to distinguish the potency of ligands.

GTM has also been applied in de novo design. For the first time in 2014, Mishima et al. [START_REF] Mishima | Development of a New De Novo Design Algorithm for Exploring Chemical Space[END_REF] used GTM for an assessment of biological activities for virtually enumerated structures. Another attempt of the compound generation with specific activity(ies) has been made with Stargate GTM [START_REF] Gaspar | Stargate {GTM}: Bridging Descriptor and Activity Spaces[END_REF]. Stargate GTM uses two manifolds: one is built in the descriptor space and another in the activities space; thus, the two spaces are bound. A specific mapping function allows to "warp" from the activities space to descriptors space, therefore identifying the values of the desirable descriptors. Once the values of the descriptors are found, it is needed to generate structures with high similarity to the detected descriptors vector, thus assuming that the generated structures would posses the searched activity. Recently GTM has been combined with auto-encoder, where the map was trained on the generated latent descriptors. Sattarov et al. [START_REF] Sattarov | De Novo Molecular Design by Combining Deep Autoencoder Recurrent Neural Networks with Generative Topographic Mapping[END_REF] used this approach to generate and analyze the binding potency of ligands of Adenosine A2a receptors. A similar approach was applied to the discovery of novel chemical reactions [START_REF] Bort | Discovery of Novel Chemical Reactions by Deep Generative Recurrent Neural Network[END_REF]. The authors have trained a sequence-to-sequence autoencoder on the USPTO [START_REF] Lowe | Chemical reactions from US patents[END_REF] reaction database. The autoencoder latent space has been visualized using GTM, the zones of the map populated by Suzuki reactions have been targeted. Many of the generated chemical reactions possessed reaction centers not present in the training set.

During this thesis, a broader exploration of GTM capabilities in chemoinformatics tasks such as virtual screening and diverse library selection has been done. Sidorov et al.

[10] have successfully built a universal map -a GTM-based multi-target classification model. In the first project of this thesis described in chapter 4, entitled "Consensus modeling using universal maps", we have applied the same model building protocol on the data extracted from ChEMBL v.23 (1.5M compounds with known activities on 618 targets).

In this project, several universal generative topographic maps have been obtained, each map being built in different descriptor spaces (hence encoding different distinct structural features). For each target-specific subset of 618 targets, the balanced accuracy (BA) has been computed; the score used to quantitatively describe the predictive performance of the map was calculated by averaging the BA over all the 618 target-specific subsets. The obtained universal maps have shown similar scores. The results are shown in chapter 4, and they answer on the following questions: i) For a virtual screening task, is it better to use one sole "best map" or several maps in consensus? ii) If the latter -how many maps should be applied?

The second project is presented in chapter 5, "In silico mining for new Bromodomain 4 inhibitors". Bromodomain 4 (BRD4) [START_REF] Sanchez | The bromodomain: from epigenome reader to druggable target[END_REF][START_REF] Liu | Drug discovery targeting bromodomain-containing protein 4[END_REF] can be considered a difficult target for virtual screening because of its flexible structure -2 α-helixes bound with two loops.

Known inhibitors of BRD4 usually form 1 H-bond with the target, the rest of the proteinligand interactions being hydrophobic [START_REF] Sanchez | The bromodomain: from epigenome reader to druggable target[END_REF]. This project was carried out in collaboration with the Enamine company [48], its goal was to find new inhibitors of BRD4. Two public sources of data have been used to form a training set -Reaxys and ChEMBL. The obtained models were used for virtual screening of 2M compounds that are available at Enamine.

Our collaborators agreed to test 3000 compounds; hence our goal was to find the top 3000 compounds that are most likely to be BRD4 inhibitors. Here, GTM was used in rather complex virtual screening funnel, including SVM models [START_REF] Chang | {LIBSVM}: A library for support vector machines[END_REF], ligand-based pharmacophores [START_REF] Wolber | LigandScout: 3-D pharmacophores derived from protein-bound ligands and their use as virtual screening filters[END_REF]50] and docking [START_REF] Hoffer | S4MPLE--Sampler For Multiple Protein--Ligand Entities: simultaneous docking of several entities[END_REF].

The third project is described in chapter 6, "In silico speciation assessment of Dynamic Combinatorial Libraries". In this project, GTM performance in diverse library selection has been compared to the classical dissimilarity-based method -MaxMin [START_REF] Holliday | Definitions of" dissimilarity" for dissimilaritybased compound selection[END_REF].

Dynamic combinatorial libraries (DCL) [START_REF] Lehn | Dynamic combinatorial chemistry and virtual combinatorial libraries[END_REF]8] are the cornerstone of the dynamic combinatorial chemistry. It relies on the reversible nature of the involved reactions between the constituents of the libraries. Briefly, a DCL is usually represented as a solution of mxn reactants that can reversibly interact one with another, leading to a formation of multiple reaction products and their thermodynamic stability dictates their distribution. An external effector (e.g., biological target) is introduced in DCL; it can lead to a global change of the distribution of the products in solution according to Le Chatelier principle if the effector will selectively bind to one (or a few) members of the DCL. It can be done in order to find the "best binder" with a given target. In this case, the reaction products have to be in almost equal proportions, because for a solution where one or a few constituents are overrepresented, the preferred interaction of a minor DCL constituent with the target may not be strong enough to overturn the equilibrium. There is software that can predict the numerical distributions of the amount of substances in the solution. However, all of them require precise thermodynamic data on equilibrium constants of the involved reactions, which complicates their usage because of the availability of the thermodynamic data for the specific case. One way to overcome this constraint is the usage of QSPR models that will predict the equilibrium constants, but again for the training of the QSPR model, the data is needed. In this project, we present the first steps toward in silico DCL modeling on the example of imine-based DCL with human Carbonic Anhydrase II as the effector. The presented workflow requires data on both equilibrium constants of the involved reactions and binding affinities with the biological target. The data on binding affinities can be extracted from public databases like ChEMBL, but there is less data on equilibrium constants of reactions. We have selected a diverse library of imines and our collaborators from prof. Lehn's laboratory synthesized and measured the equilibrium constants for over 250 reactions corresponding to the selected pool of imines, which became the training set in this study.

Methods

In this section of the thesis, several aspects of the applied methods will be discussed.

First of all, the basics of QSAR methodology are described, which include a brief description of the QSAR paradigm, some words on different types of descriptors and popular machine learning methods used in model building. It will be followed by a description Support Vector Machine (SVM) and Generative Topographic Mapping (GTM) used as machine-learning methods in the presented projects.

QSAR /QSPR methodology

Quantitative structure-activity relationship (QSAR)/Quantitative structure-property relationship (QSPR) modeling is one of the cornerstones of the chemoinformatics tools regularly used in many fields such as medicinal chemistry and material science [START_REF] Cherkasov | QSAR modeling: where have you been? Where are you going to?[END_REF]. The principle of the modeling is to find a mathematical function that relates a chemical structure to the studied property (such as logP) or activity (for a given biological target). Such modeling implies that the structural information of the molecule is encoded in numerical form -in a vector of molecular descriptors [START_REF] Todeschini | Substructural fragments: an universal language to encode reactions, molecular and supramolecular structures[END_REF], and the values of these molecular descriptors are used to define the position of the compounds in the chemical space. All of the above said could be expressed by the following equation:

= ( )

To model any possible property/activity, one needs a way to encode all (or at least most of) the essential structural information and a big enough dataset of compounds with known activity. Last but not least, a machine learning algorithm that will eventually build a predictive model.

When it comes to the calculation of molecular descriptors, there are several ways to do it. For instance, 1D molecular descriptors are directly obtained from the chemical formula of the compound. These descriptors are the most straightforward descriptors related to the molecule's "fundamental properties," such as the number of atoms and (central atoms with their environment) or triplets that encode the compound's atoms and/or bond types. In addition to this, the fragments may be colored, adding some additional information: pharmacophoric types of atoms, formal charges force-field atom types, etc.

The dataset that will be used for model training [57] is called the training set. The models can be divided into two groups basing on the modeled property or activity:

regression model (when the modeled/predicted property is a numerical value) and classification model (when the modeled/predicted property is categorical). The machine learning algorithms that require a known (experimental) property are called supervised learning algorithms. A non-exhaustive list of these methods includes Multilinear Regression (MLR) [START_REF] Marill | Advanced statistics: linear regression, part II: multiple linear regression[END_REF], Random Forest (RF) [START_REF] Liaw | Classification and regression by randomForest[END_REF], Support Vector Machine (SVM) [START_REF] Cortes | Support-vector networks[END_REF] and Artificial Neural Network (ANN) [START_REF] Hassoun | Fundamentals of artificial neural networks[END_REF]. Supervised modeling occurs when the modeled property/activity is known for all the entries in the training set, and the training per se consists of fitting a function in order to minimize the prediction error. Another type of machine-learning algorithms is unsupervised learning, where the compound's property is not used. Usually, the unsupervised methods are used in order to extract some "inherited" structural information basing on data distribution, to describe and interpret the data and, in some cases, to visualize it. This class of methods includes clustering (hierarchical clustering [START_REF] Johnson | Hierarchical clustering schemes[END_REF], k-means [START_REF] Lloyd | Least squares quantization in PCM[END_REF]) and dimensionality reduction (Principal Components Analysis (PCA) [START_REF] Hotelling | Analysis of a complex of statistical variables into principal components[END_REF], Self-Organizing Maps (SOM) [START_REF] Kohonen | The self-organizing map[END_REF], Generative Topographic Maps (GTM) [START_REF] Bishop | GTM: The generative topographic mapping[END_REF]).

An important parameter of any model is its quality [START_REF] Leach | An introduction to chemoinformatics[END_REF]. Usually, for regression models, Root Mean Squared Error (RMSE) and determination coefficient (R 2 ) are used.

RMSE is calculated according to the formula below, where N is the number of compounds, yexp,i and ypred,i are experimental and predicted property values of ith molecule respectively:

= ∑ , - , R 2
estimates the correspondence of experimental and predicted values. The maximal value of R 2 is 1, which corresponds to an ideal fit of the model, i.e., all the predicted values are equal to experimental ones. Lower values of R 2 correspond to a worse model and "acceptable" value of R² > 0.5 [START_REF] Golbraikh | Beware of q2![END_REF]. Determination coefficient is calculated by the formula below, where yexp,i and ypred,i are experimental and predicted property values of ith molecule respectively and <yexp> is the average property value of the dataset:

= 1 - ∑ , - , ∑ , -〈 〉
When it comes to classification models, the model's quality is related to the number of compounds with the correctly assigned category. Usually, the classification tasks are reduced to binary classification, where the compounds are split into 2 classes (active/inactive). In this case, to evaluate the model's quality, a confusion matrix is used. It represents a table where the predicted class of the compounds is matched with their actual value. The cells of the table contain the number of compounds that have had a correct ("True") or an incorrect ("False") class assignment that is denoted as "Positive" and "Negative". Using the confusion matrix, one can calculate the balanced accuracy (BA), which is a numerical characteristic of a classification model. BA takes the rate of correct predictions of both classes in equal proportions, and it takes values from 1 (ideal case) to 0.5 (random predictions). BA is used for the datasets where the predomination of one class over another is observed:

Actual

= 1 2 + + +
Another essential factor that is intrinsically related to model quality is its applicability domain (AD). The machine-learning methods perform well in interpolating tasks because the model is trained on a limited set of compounds. The model provides reliable predictions for the compounds that are similar to the ones from the training set. In every project of this thesis, ISIDA descriptors were used; fragment control was generally used as AD [START_REF] Varnek | ISIDA-Platform for virtual screening based on fragment and pharmacophoric descriptors[END_REF]. This approach considers any compound to be out of the AD if it has at least one descriptor (i.e., substructural molecular fragment) that was not present within the training set.

One cannot tell how well a model will perform in a "real-life scenario", nonetheless basing on its performance during cross-validation (CV) [START_REF] Browne | Cross-validation methods[END_REF], one can estimate the model's

Support Vector Machine

Support Vector Machine (SVM) is a popular supervised machine learning method that can be used in classification and regression tasks. The method has been developed and published by Vapnik in 1995 [START_REF] Cortes | Support-vector networks[END_REF] following the idea to find a hypersurface that separates two classes of objects with the "gap" (called margin) as wide as possible. New objects that will be mapped to the same space will be assigned to one of the classes according to the side of the surface where they fall on. Additionally, SVM can non-linearly map the input vectors into higher dimensional feature space using a kernel function (the so-called kernel trick) and then to linearly separate them in this new feature space. However, perfect separation of classes is not possible, since "an ideal SVM" should produce a hyperplane able to completely separate the objects of different classes. Such hyperplane may result in an overfit model; therefore, new objects might be wrongly classified. The SVM algorithm is maximizing the margin, and in the meantime, it minimizes the misclassifications using the slack variable ξi. The goal of the algorithm becomes to maintain the slack variable at a minimal value while maximizing the margin; therefore the constraint and objective function becomes

( • + ) 1 -, ∀ 0 min ‖ ‖ + ∑ ,

with C being trade-off margin

In 1996 Vapnik proposed a support vector regression (SVR) [START_REF] Drucker | Support vector regression machines[END_REF] as the development of the SVM method for the prediction of a continuous variable. SVR maintains the main feature that characterizes the original algorithm -the maximal margin. In cases of regression task ( = • + ), a margin of error tolerance ε is set in approximation to the hyperplane and the algorithm is minimizing:

1 2 ‖ ‖ + ( + * )
Where C is the cost that defines the penalty for objects whose predicted value deviates from the experimental value for more than ε. The constraints become:

-( • + ) ≤ ε + ( • + ) - ≤ ε + * * , 0

Generative Topographic Mapping

Generative Topographic Mapping (GTM) has been presented by Bishop in 1998 [START_REF] Bishop | GTM: The generative topographic mapping[END_REF] as a method of data visualization. GTM is a probabilistic extension of Self-Organizing Maps (SOM) [START_REF] Kohonen | The self-organizing map[END_REF], but unlike SOM, it considers the likelihood of the training data as the objective function. Moreover, in GTM, a single object is not associated with one particular node (and its neighbors), but it is associated with the probability distribution over the entire latent space.

GTM finds a representation of the data distribution in the initial D-dimensional data space on an L-dimensional hypersurface called a manifold. Although the dimensionality of the manifold is usually L=2, L can take any natural value from 1 to D. To map the objects from the initial space to the latent space, a mapping function y(x, W) given as a grid of M Gaussian activation functions (radial basis functions, or RBF) is applied:

( , ) = exp ‖ -‖ 2
Where M is the number of RBFs, D is the initial space dimensionality, W is the matrix ( × ) of weights connecting the initial data space and RBF grid, xm is the center of the m-th RBF; d takes values from 1 to D, M and σ being the parameters of the method.

Since the function y(x, W) is smooth (and therefore continuous), the so-called neighborhood behavior is observed -the objects that are close in the initial space remain neighbors in the latent space. Every node of the grid is associated with the center of a normal distribution function with inverse variance β, that corresponds to the sampling of the random variable t with the following probability density function:

( | , ) = 1 exp -2 ‖ -( , )‖
Where K is the number of nodes, xk the coordinate of the k-th grid node in the latent space, ( , ) -the coordinates to which it has been mapped in the initial data space and t covers the whole data space representing any object. The logarithm of the probability with which the data could be generated is called log-likelihood, and it is denoted as LLh. The higher the value of LLh, the better the manifold represents the data. LLh is used as the maximized function in the Expectation-Maximization (EM) algorithm. LLh is a function of two parameters, W and β:

LLh( , ) = LLh = ln 1 exp -2 ‖ -( , )‖
tn is the position of the n-th object in the initial data space. Every object has a non-zero probability of being mapped into any node of the grid. This probability is called responsibility, and it is calculated using Bayes's theorem:

= ( | ) = exp - 2 ‖ -( , )‖ ∑ exp - 2 ‖ -( , )‖
For every object, the responsibility is normalized over the grid of nodes; therefore, the sum of responsibilities for a given object is 1. This vector is used for visualization and modeling purposes.

Standard GTM (sGTM) [START_REF] Gaspar | Generative Topographic Mapping Approach to Chemical Space Analysis[END_REF] algorithm is computationally expensive in the case when the number of the compounds and descriptors is relatively high since it takes too long to calculate the Euclidean distances for all the pairs "object-node". For a small dataset of 1000 compounds having 500 descriptors, the computing time of all the distances between all the compounds and 900 nodes is around 4.5s on a single CPU (Intel Core i7-6700HQ) [START_REF] Lin | Parallel Generative Topographic Mapping: an Efficient Approach for Big Data Handling[END_REF].

However, the computational time rises to 135s for the dataset of just 30k compounds. Since this procedure is done at each iteration of the EM algorithm, it renders the sGTM algorithm rather slow. Moreover, when it comes to extensive collections of data, a memory problem might appear. To overcome the constraints of sGTM, an incremental variation of GTM (iGTM) [START_REF] Gaspar | Chemical data visualization and analysis with incremental generative topographic mapping: big data challenge[END_REF] has been proposed. In this case, the initial dataset, instead of being processed as a whole, is divided into many blocks, then each block is processed consecutively. The manifold will be trained on one single block at a time, which accelerates the procedure. In the current work, iGTM has been used almost exclusively.

Recently a new approach in GTM methodology has been proposed -parallel GTM (pGTM) [START_REF] Lin | Parallel Generative Topographic Mapping: an Efficient Approach for Big Data Handling[END_REF]. The main idea of pGTM is to extend iGTM over multiple CPUs. It is done by manifold initialization over all dataset, then this dataset is split into blocks, and the initialized manifold is fit on each block independently and simultaneously. Thus, each block provides an intermediate manifold that has been fitted on a given part of the data. In the end, an averaging of the W and β over all blocks is done.

GTM as a visualization method and modeling tool

For the visualization and modeling of the data, GTM uses the concept of landscape.

For every compound, GTM generates a vector of normalized responsibilities that can be treated in the same way as molecular descriptors. The number of these descriptors will be equal to the number of the nodes used in the GTM grid. The landscape is obtained by adding a specific property of interest for the given dataset as a third dimension of the 2D map. Three types of landscapes can be defined: class landscape [START_REF] Gaspar | Generative topographic mapping-based classification models and their applicability domain: application to the Biopharmaceutics Drug Disposition Classification System (BDDCS)[END_REF], property landscape [START_REF] Gaspar | GTM-Based QSAR Models and Their Applicability Domains[END_REF] and density landscape [START_REF] Lin | Parallel Generative Topographic Mapping: an Efficient Approach for Big Data Handling[END_REF].

The Where Ntot is the total number of training items, and rkn is the responsibilities of the members of class ci in the node k. To predict a class for a new compound q, the following equation is used:

( | ) = ∑ ( ) ( | ) = ( | ) × ( ) ∑ × ( ) =
= ( | ) ×
To visualize a class landscape, the normalized probability of class ci is used as a 3 rd axis (color code). The population of the nodes is taken into account by the addition of the transparency to the used colors.

GTM-based regression models are relying on the property landscapes, which represent the distribution of a property over the latent space. The definition of property landscape is done by using a list of property values of compounds that correspond to a particular node:

= ∑ × ∑
Where pn is the property value for the compound n, and pk is the mean property value for the node k. The visualization of the property landscape is done by using the pk values as the 3 rd axis and interpreting them as color code. In contrast, transparency is used (same as in the class landscapes) to take into account the nodes population.

The property of a new compound q is predicted similarly to the class prediction:

= ×

The density landscape could be viewed as a "subtype" of property landscape, where pk represents the sum of all compound responsibilities in the node k. This landscape is usually applied for the analysis of the data distribution when there is no particular property to use (see chapter 6.1.2.2 Cell-based diverse-library selection using GTM) or when the landscape transparency is not easily readable.

Applicability domain of GTM-based QSAR models

While using GTM for modeling, one can use classical AD approaches, but GTM is offering several definitions of AD. Several different GTM-based AD definitions have been reported [START_REF] Gaspar | GTM-Based QSAR Models and Their Applicability Domains[END_REF]: likelihood-based, density-based and class-dependent density.

When applying likelihood-based AD, a compound is considered to be out of the GTM AD if its position is too far away from the manifold in the initial data space. To apply this AD concept, the LLh cutoff is determined by sorting the training compounds accordingly to their LLh from the higher to lower LLh value. The cutoff is set at n% of compounds (usually 5%) having the smallest LLh; thus, the LLh cutoff is taken as the highest LLh out of this "bottom" n%.

The density-based AD discards the nodes on the GTM landscape, where the cumulative responsibility is below a certain threshold. This AD allows using only populated zones to make the predictions. This concept is the easiest to visualize and interpret since all the testing compounds that have been projected on the white zones of the map are considered to be out of AD.

The class-dependent density AD is similar to the density-based AD. The difference is that the density of the winning class cbest in the node is checked, which has the highest conditional node probability P(xk|cbest). The ratio for this predominance is a user-defined variable. This AD concept is especially useful for classification models.

Universal GTM

Universal Generative Topographic Maps (uGTM) are GTM-based classification models that have been introduced by Sidorov et al. [START_REF] Sidorov | Mappability of drug-like space: towards a polypharmacologically competent map of drug-relevant compounds[END_REF]. In this work, the authors aimed to cover large chemical space defined by the ChEMBL database of version 20 using a single map. The main difference between uGTM and "local" GTM (GTM explicitly built for a dataset of compounds manifesting a particular activity) is the data used for manifold fitting.

Usually, when one is applying GTM, the compounds of the training set share a common activity/property; thus, the obtained map will describe its landscape. The data used for uGTM contains more than 1.2M ligands with known activity for more than 400 biological targets. The descriptors space and the GTM parameters were selected using the Genetic algorithm [START_REF] Gaspar | Generative Topographic Mapping Approach to Chemical Space Analysis[END_REF][START_REF] Horvath | An evolutionary optimizer of libsvm models[END_REF] described in chapter SVM/GTM parameters tuning. The results showed that the uGTM approach could efficiently cover a broad range of chemotypes. The best map selected by GA was cross-validated on 410 ChEMBL targets, showing that approximately 80% of the targets were predicted with the mean Balanced Accuracy of 0.7.

SVM/GTM parameters tuning

The performance of the machine learning methods is parameter dependent. For instance, the SVM/SVR performance depends on the type of kernel and the regularization coefficient C, while the GTM has four parameters: number of nodes k, number of RBFs m, regularization coefficient l, and RBF's width w. Besides these parameters, an "optimal" descriptor space is also needed to be found. To tune all these parameters, the genetic algorithm (GA) has been used in all of the projects. GA is a stochastic approach that allows achieving "the maximal" model performance while trying a variety of combinations of parameters from a pre-defined range.

The algorithm's details have already been described in several publications [START_REF] Gaspar | Generative Topographic Mapping Approach to Chemical Space Analysis[END_REF][START_REF] Horvath | An evolutionary optimizer of libsvm models[END_REF].

Shortly, GA generates a set of chromosomes; each chromosome is presented by a vector of model's parameter values, as well as some meta-parameters like descriptor space. Each attempt (chromosome) is validated using n-fold cross-validation repeated m times, and a fitness score (FSc) is associated with the attempt. The FSc is related to the model's success with a given chromosome; for classification tasks, FSc is related to BA and for regression tasks to R 2 . Higher scored chromosomes will be allowed to generate "children" using crossovers and mutations, which might result in potentially better FSc. The GA stops in two cases: either no FSc improvement has been observed during the last two generations, or the maximal number of attempts has been achieved.

In the case of SVR-or GTM-based regression model, the FSc is defined by a crossvalidated determination coefficient. For each repetition of the n-fold cross-validation the R is computed. By default, the algorithm is doing 3-fold cross-validation repeated 12 times. Then, the mean value of R (<R 2 >) as well as its standard deviation σ is calculated.

The FSc is defined as:

= 〈 〉 -2 × σ
When the GA is run for the optimization of the classification model, the FSc is calculated as follows:

= 〈 〉 -2 × σ
Where <BA> is the mean value of the cross-validated BA of each repetition and σ being its standard deviation.

4 Consensus modeling using universal maps

Introduction

Virtual Screening (VS) [START_REF] Reddy | Virtual screening in drug discovery-a computational perspective[END_REF] is a technique applied in drug discovery to search libraries of molecules in order to identify the compounds with the property/activity of interest using knowledge retrieved from the existing data. Usually, the so-called VS funnel has several layers of applied methods differentiating them in terms of accuracy/speed ratio.

For instance, the methods having low accuracy but high computational speed (like filters or similarity search) will usually be applied in the first place in order to eliminate the compounds that are less likely to be active. On the other hand, methods providing high accuracy with the cost of slower computational speed (like docking) will be applied at "the bottom of the funnel" on a more restricted set of compounds since they are more likely to be active.

GTM has proven to be able to produce target/property-specific models having comparable performance to other machine learning methods like SVM and RF. However, in contrast to these methods, the manifold fitting is an unsupervised process. Therefore, with GTM, one manifold can fit any database containing thousands or even millions of compounds. These compounds may have various activities/properties; hence with only one fitted manifold, one have access to all the landscapes representing the present compounds' activities/properties. This has been applied and tested in the work of Sidorov et al. [START_REF] Sidorov | Mappability of drug-like space: towards a polypharmacologically competent map of drug-relevant compounds[END_REF]., where for the first time, the concept of universal GTM (uGTM) has been introduced.

Here, the same protocol [START_REF] Sidorov | Mappability of drug-like space: towards a polypharmacologically competent map of drug-relevant compounds[END_REF] has been applied to data extraction (ChEMBL 23) and standardization, as well as uGTM generation. A total number of 1.5M compounds with known activities on 618 targets have been extracted. The subsets of the ligands of 236 targets, including GPCRs, kinases, nuclear receptors etc., have been used for 3-fold crossvalidation, and 382 target-specific subsets have been used exclusively as external-validation sets. Directory of Useful Decoys (DUD) [START_REF] Huang | Benchmarking Sets for Molecular Docking[END_REF] has been used as a genuinely externalvalidation set. The DUD dataset has been standardized in the same way as previously extracted compounds from ChEMBL. The standardization has been followed by removal from the DUD dataset of the compounds that are already present in ChEMBL in order to create orthogonal external data sets. Most targets had a complete overlap of active compounds when they were simultaneously present in ChEMBL and DUD. In these cases, corresponding target-specific sets have been discarded, however, in nine cases the DUD database contained sufficiently numerous original actives, thus leading to 9 target-specific subsets. Eight uGTMs have been selected being ranked as "best maps" by genetic algorithm, with the scoring function being the average BA for all target-specific landscapes on 3-fold cross-validation. Each uGTM is based on different ISIDA descriptor space, encoding distinct structural features. Although the average BA of eight maps is roughly equivalent, the maps were showing different BAs on target-specific subsets. For instance, an individual map was showing high BA value for a given target-specific subset while having a relatively low value on another target-specific subset; moreover, another map could show the opposite behavior. This induced an in-depth analysis of 8 selected uGTMs.

Performance evaluation of universal maps

The uGTMs performance was evaluated using three scores: i) BA in 3-fold CV (using ChEMBL compounds) and in VS (using DUD); ii) Receiver Operating Characteristic Area Under Curve (ROC AUC) in VS; iii) Enrichment Factor (EF) in VS. BA has been mainly used during cross-validation. BA serves to assess the ability of landscapes to predict the correct activity class of candidates not used for landscape construction, i.e., both in "internal" cross-validation and "external" VS. Note that reported BA scores for individual maps -both in the CV and in VS applications -are always calculated on the entire targetspecific sets. It includes all DUD compounds, even those projected onto empty map zones; hence these molecules are out of AD, and they are considered, by default, inactive.

However, ROC AUC is a more natural VS evaluation criterion than BA, since, in VS, the critical element is the relative ranking of candidates -a significant prioritization of the active compounds with respect to the inactive. The ranking was performed according to the GTM landscape-predicted probability of each compound to be active. The compounds falling outside the applicability domain were assigned zero probability of activity; thus, they were placed at the bottom of the ranking list. EF for the top 100 ranked molecules was calculated according to the equation below:

Conclusions

In this work, the predictive performance of eight newly constructed uGTM models in a "strict" 3-fold CV and VS of nine target-specific subsets of compounds extracted from DUD has been assessed. It has been shown that these maps can provide a relatively good separation (BACV >0.6) of active and inactive for the majority of 618 ChEMBL targetspecific subsets, irrespective of whether these subsets have been used in model training or not. It has been found out that any individual map could not achieve consistently accurate predictions for each target-specific subset. However, it has been proven that these maps, which were each built on a different descriptor space, are highly complementary -the target-specific series of compounds that are being predicted poorly by one uGTM will be much better predicted by another. For 617 out of 618 activity classes, at least one uGTM provides a highly discriminatory activity landscape.

It was observed that there is no correlation between performance in the CV and external predictive power of individual activity landscapes. A solution has been found -a consensus approach. The most important advantages of this approach are 1) 100% data coverage in most of the cases; 2) a significant increase in EF for the 100 top-ranked compounds; 3) high performance of the consensus model compared to individual models based on ROC AUC. Last but not least, seven uGTMs have been proven to be sufficient to provide complementary views of biologically relevant chemical space that resulted in the enhancement of the performance in VS. 5 In silico mining for new Bromodomain inhibitors

CHEMBL3202 CHEMBL3351 CHEMBL3535 CHEMBL3721 CHEMBL5524 CHEMBL321 CHEMBL3356 CHEMBL3864 CHEMBL3729 CHEMBL5543 CHEMBL3227 CHEMBL3357 CHEMBL3869 CHEMBL3746 CHEMBL5545 CHEMBL3230 CHEMBL3359 CHEMBL3880 CHEMBL3759 CHEMBL5568 CHEMBL3385 CHEMBL3589 CHEMBL3764 CHEMBL3886 CHEMBL6003 CHEMBL3397 CHEMBL3590 CHEMBL3772 CHEMBL3890 CHEMBL6007 CHEMBL3399910 CHEMBL3616 CHEMBL3776 CHEMBL3891 CHEMBL6154 CHEMBL340 CHEMBL3622 CHEMBL3778 CHEMBL3892 CHEMBL4895 CHEMBL3401 CHEMBL3629 CHEMBL3785 CHEMBL3898 CHEMBL4896 CHEMBL3426 CHEMBL3636 CHEMBL3788 CHEMBL3902 CHEMBL4897 CHEMBL3437 CHEMBL3650 CHEMBL3795 CHEMBL3905 CHEMBL4898 CHEMBL3438 CHEMBL3663 CHEMBL3807 CHEMBL3906 CHEMBL4899 CHEMBL3468 CHEMBL3683 CHEMBL3816 CHEMBL3911 CHEMBL4444 CHEMBL3474 CHEMBL3687 CHEMBL3819 CHEMBL3913 CHEMBL4461 CHEMBL3475 CHEMBL3691 CHEMBL3820 CHEMBL3920 CHEMBL4462 CHEMBL3476 CHEMBL3961 CHEMBL3829 CHEMBL3922 CHEMBL4465 CHEMBL3510 CHEMBL3965 CHEMBL3831 CHEMBL3935 CHEMBL4478 CHEMBL3514 CHEMBL3969 CHEMBL3835 CHEMBL3959 CHEMBL4481 CHEMBL3836 CHEMBL3972 CHEMBL4051 CHEMBL4203 CHEMBL4482 CHEMBL3837 CHEMBL3973 CHEMBL4068 CHEMBL4204 CHEMBL4501 CHEMBL3861 CHEMBL3974 CHEMBL4071 CHEMBL4223 CHEMBL4506 CHEMBL3863 CHEMBL3975 CHEMBL4072 CHEMBL4224 CHEMBL4801 CHEMBL3572 CHEMBL3976 CHEMBL4073 CHEMBL4225 CHEMBL4803 CHEMBL3582 CHEMBL3979 CHEMBL4079 CHEMBL4227 CHEMBL4804 CHEMBL3587 CHEMBL3982 CHEMBL4080 CHEMBL4234 CHEMBL4816 CHEMBL3983 CHEMBL4081 CHEMBL4237 CHEMBL4422 CHEMBL4581 CHEMBL3991 CHEMBL4093 CHEMBL4247 CHEMBL4426 CHEMBL4599 CHEMBL4005 CHEMBL4101 CHEMBL4261 CHEMBL4427 CHEMBL4600 CHEMBL4015 CHEMBL4123 CHEMBL4270 CHEMBL4439 CHEMBL5261 CHEMBL4016 CHEMBL4128 CHEMBL4273 CHEMBL4441 CHEMBL5282 CHEMBL4018 CHEMBL4142 CHEMBL4282 CHEMBL4714 CHEMBL5285 CHEMBL4026 CHEMBL4145 CHEMBL4296 CHEMBL4718 CHEMBL5314 CHEMBL4029 CHEMBL4147 CHEMBL4302 CHEMBL4722 CHEMBL5330 CHEMBL4036 CHEMBL4158 CHEMBL4303 CHEMBL4761 CHEMBL5331 CHEMBL4040 CHEMBL4176 CHEMBL4306 CHEMBL4766 CHEMBL6164 CHEMBL4045 CHEMBL4179 CHEMBL4315 CHEMBL4608 CHEMBL6166 CHEMBL4374 CHEMBL4191 CHEMBL4338 CHEMBL4617 CHEMBL6175 CHEMBL4375 CHEMBL4198 CHEMBL4361 CHEMBL4618 CHEMBL4698 CHEMBL4376 CHEMBL4202 CHEMBL4367 CHEMBL4625 CHEMBL4699 CHEMBL4393 CHEMBL4508 CHEMBL4662 CHEMBL4630 CHEMBL4852 CHEMBL4394 CHEMBL4516 CHEMBL4674 CHEMBL4576 CHEMBL4829 CHEMBL4398 CHEMBL4523 CHEMBL4681 CHEMBL4578 CHEMBL4835 CHEMBL4408 CHEMBL4525 CHEMBL4683 CHEMBL4708 CHEMBL4601 CHEMBL4822 CHEMBL4575 CHEMBL4685

Introduction

The goal of this project was to carry out a virtual screening (VS) of a dataset of 2M compounds provided by Enamine company [48], in order to find new inhibitors of Bromodomain 4 (BRD4) [START_REF] Liu | Drug discovery targeting bromodomain-containing protein 4[END_REF]. The provided dataset of 2M compounds corresponds to the compounds that are physically available in stocks, or that could be easily synthesized. It has been agreed to provide our Enamine collaborators with a dataset of 3000 compounds that would be tested. To get this dataset of 3000 compounds, a VS screening protocol has been developed, which included a consensus application of GTM and SVM classification models as well as ligand-based pharmacophore models.

Bromodomain 4

Biological role

Lysine acetylation of histone proteins is a fundamental post-translational modification regulating chromatin structure, and it plays a significant role in gene transcription [START_REF] Kouzarides | Chromatin modifications and their function[END_REF].

Readers of post-translational modifications are structurally diverse proteins than contain one or more effector modules that recognize covalent modifications of proteins and DNA.

The recognition of acetylation of lysine residues is primarily initiated by bromodomains [START_REF] Sanchez | The bromodomain: from epigenome reader to druggable target[END_REF]. Bromodomains are involved in the regulation of transcriptional programmers. They have been identified in oncogenic rearrangements [START_REF] Jain | Bromodomain histone readers and cancer[END_REF] that lead to highly oncogenic fusion proteins, which play a crucial role in the development of several aggressive types of cancer [START_REF] Sanchez | The bromodomain: from epigenome reader to druggable target[END_REF][START_REF] Stathis | BET proteins as targets for anticancer treatment[END_REF] (like NUT carcinoma, leukemia and lymphoma [START_REF] Marcotte | Functional genomic landscape of human breast cancer drivers, vulnerabilities, and resistance[END_REF][START_REF] French | NSD3--NUT fusion oncoprotein in NUT midline carcinoma: implications for a novel oncogenic mechanism[END_REF]).

Recently it has been shown what role is playing BRD4 in NUT carcinoma. NUT carcinoma is a very aggressive and rare form of undifferentiated squamous-cells carcinoma.

It is considered one of the most lethal solid tumors, which typically is non-responsive to chemotherapy or radiotherapy and an overall survival spanning from 6 to 9 months [START_REF] French | BRD--NUT oncoproteins: a family of closely related nuclear proteins that block epithelial differentiation and maintain the growth of carcinoma cells[END_REF].

This disease is genetically defined by chromosomal rearrangements involving the NUT gene fused to the BRD4 [START_REF] French | BRD--NUT oncoproteins: a family of closely related nuclear proteins that block epithelial differentiation and maintain the growth of carcinoma cells[END_REF]. This creates BRD4-NUT oncogene that is considered to be a main pathogenetic driver of cellular transformation. It has been found that the interception of the BRD4-NUT fusion gene results in the slowing of the differentiation and growth of NUT carcinoma cells [START_REF] French | NSD3--NUT fusion oncoprotein in NUT midline carcinoma: implications for a novel oncogenic mechanism[END_REF][START_REF] Filippakopoulos | Selective inhibition of BET bromodomains[END_REF].

BRD4 as a therapeutic target

Bromodomain modules share a conserved fold that comprises a left-handed bundle of four α-helices (named αZ, αA, αB and αC) [START_REF] Sanchez | The role of human bromodomains in chromatin biology and gene transcription[END_REF] that are linked by diverse loop regions of variable charge and length (known as ZA and BC loops) which surround a central acetylated lysine binding site (Figure 5.1). BRD4 can be considered a difficult target for virtual screening because of its flexible structure. Known inhibitors of BRD4 usually form 1 hydrogen bonds with the protein, and the rest of the protein-ligand interactions being hydrophobic [START_REF] Sanchez | The bromodomain: from epigenome reader to druggable target[END_REF].

Methods

The virtual screening protocol (Figure 5

Pharmacophore models

IUPAC's definition of the pharmacophore model is "an ensemble of steric and electronic features that is necessary to ensure the optimal supramolecular interactions with a specific biological target structure and to trigger (or to block) its biological response" [START_REF] Ganellin | Glossary of terms used in medicinal chemistry[END_REF].

A pharmacophore does not represent a real molecule or a real association of functional groups, but a purely abstract concept that accounts for the common molecular interaction capacities of a group of compounds towards their target structure. The pharmacophore can be considered as the largest common denominator shared by a set of active molecules. The pharmacophore features are H-bond acceptors and donors, charged or ionizable groups, hydrophobic groups and aromatic rings.

The spatial relationships between the features in a 3D pharmacophore model can be specified as distances or distance ranges or by defining the (xyz) locations of the features together with some distance tolerance (typically as a spherical tolerance region). There are different possibilities to derive pharmacophore models: based on the three-dimensional structure of a ligand-protein complex (structure-based modeling) or based on the structural information of active compounds only (ligand-based modeling). In this project, ligandbased pharmacophore models have been developed.

Generally, a database is built in such a way that the molecules that it contains are usually (or at least it is expected) represented by a set of conformers that supposedly include the bioactive geometry adopted during the interaction with the target protein. All conformers of used compounds are superimposed, and the associated common pharmacophore features are generated. Then, it is up to the user to define the number and the types of needed pharmacophores that will form the model. Depending on the selectivity of the pharmacophore model, such a virtual screening of chemical databases consisting of millions of small molecules can result in tens to thousands of hits. For the compounds ranking and to model's quality determination, the matching between the pharmacophore model and each molecule of the virtual screening hit list, a score is calculated.

LigandScout [START_REF] Wolber | LigandScout: 3-D pharmacophores derived from protein-bound ligands and their use as virtual screening filters[END_REF]50] was used in the current work. The main feature of LigandScout is the fast alignment algorithm due to the efficiency of the implementation and the advanced geometric similarity measure for the chemical features. In this algorithm, the first step concerns the generation of the 3D pharmacophore features (Figure 5.3) identified for each database conformer.

Then, the algorithm creates for each feature type a set of inter-feature distances. The distance sets created for the pharmacophore model, and the conformer pharmacophore features are then compared in a pairwise manner. In order to perform a pair assignment, the so-called Hungarian matching algorithm is executed. Finally, the feature distances between model and conformer are minimized using Kabsch alignment algorithm. For estimation of alignment quality, the pharmacophore fit score function has been used [START_REF] Wolber | LigandScout: 3-D pharmacophores derived from protein-bound ligands and their use as virtual screening filters[END_REF]:

= 9 -3 × ( , 3)

= × +

Where SRMS is the matched feature pair Root Mean Squared Deviation (RMSD) score in range [0,9]; RMSFP is the RMSD of the matched feature pair distances; SFCR is the feature count/RMS distance score; c is a weighting factor for the number of matched feature pairs, and NMFP is the number of geometrically matched feature pairs.

Docking

The docking part of the project was done using S4MPLE (Sampling For Multiple Protein-Ligand Entities) software [START_REF] Hoffer | S4MPLE--Sampler For Multiple Protein--Ligand Entities: simultaneous docking of several entities[END_REF]. It is based on a hybrid genetic algorithm, which allows the simulation of one molecule (conformer generation) or many molecules (docking).

Energy calculations are done using the AMBER force field [START_REF] Wang | Development and testing of a general amber force field[END_REF] for biological macromolecules and its generalized version -GAFF for ligands. The ability of S4MPLE to indiscriminately handle inter-and intramolecular degrees of freedom is achieved through the appropriate design of torsional angles, rotational and translational degrees of freedom.

In S4MPLE, a genetic operator works on some randomly chosen covalently connected (or not) molecular substructure. If the structure is covalently connected, then the operator will affect the structure in such a way that bond length and valence angles will not be changed.

If the structure is not covalently connected, then it might be, for example, one of the ligands competing for the binding site. In that case, the guidance role of missing covalent bond will be taken by a potentially favorable contact axis, which is randomly chosen as a pair of atoms (one atom belonging to the external partner and another to structure itself), that should be brought together in order to form a hydrogen bond or a hydrophobic interaction.

The following steps make the preparation of the active site: protein atoms have to be fixed by enumerating their sequence numbers. A predefined cutoff for non-bonded interactions was established on 12Å. Protein atoms that are too far from the active site in order to ever come within 12Å to any ligand atom would merely slow down calculations by requesting the regular update of their distances to ligand atoms. Therefore, docking was not run on the entire protein, but on the selection of relevant residues that have at least one atom at less than 10Å from any of the co-crystallized ligand, herewith used to define the active site region. Moreover, S4MPLE requires the user-specified input of ''hot spots'' -key solvent-accessible atoms, chosen preferentially at the bottom of the site cavity, which will be used for random prepositioning of the ligand into the active site. These may, but do not have to, include site atoms seen to make contacts to the co-crystallized PDB ligand.

Ligands, initially provided as standardized SMILES, preprocessed by the standardization tool of the Strasbourg virtual screening web server, undergo an automated conversion, using an in-house tool developed based on the ChemAxon API, to a fully protonated initial 3D structure. The tool relies on the tautomer and, respectively, pKa plugin to generate the most probable microspecies of the expected main tautomeric form.

Users might request several tautomeric or protonation states to be generated, and each to be docked as an independent candidate -but the option was not used here. Explicit hydrogens are assigned, and the conformer plugin then generates a single conformer. Eventually, the charge plugin is used to assign Gasteiger charges [START_REF] Gasteiger | Iterative partial equalization of orbital electronegativity-a rapid access to atomic charges[END_REF] to this structure. Last but not least, the tool detects flexible rings and proposes, for each, the single bond to be formally ''broken'' in order to enable intra-cyclic torsional axes to be driven by S4MPLE.

S4MPLE docking begins by extracting all the data of a given ligand into a dedicated directory, then running a 200-generation evolutionary conformational search with S4MPLE, on the free ligand, at default settings. Next, active site data are added to this directory, and a brief fist simulation is run in order to calibrate the optimal cutoff for the interaction fingerprint dissimilarity value (minfpdiff), representing the threshold at which two conformers are considered as redundant, and thus pruned during the evolutionary process.

The proper management of population diversity has been noted to be of paramount importance for ensuring the convergence/reproducibility of evolutionary simulations. As ligands vary in sizes, so does their interaction fingerprint, making it challenging to come up with a universally applicable threshold value -hence, the need to calibrate it for each system. The population initialization procedure, regularly serving as the first step for the evolutionary simulation, is called repeatedly (10 times). After each call, the interaction fingerprints of the randomly generated population members are compared to each other, generating the complete Hamming distance matrix for all pairs of conformers in the population. The lowest, mean and maximal Hamming distances for each population are memorized. The minfpdiff threshold is defined as 90% of the average of the ten lowest intra-population Hamming distances.

Eventually, the main docking simulation is started with the above-determined minfpdiff value as a population diversity control parameter. Top poses are generated and stored together with their energy values 〈 @ 〉. The docking index ∆E $for the current ligand can be directly estimated as 〈 @ 〉 -〈 〉. After completion of docking calculations for all ligands, these can be ordered by increasing ∆E, and the final ROC curve can be generated in order to determine the area under it, as the final benchmarking criterion. The variation of the ROC AUC as a function of the performed number of generations may be informative about the minimal required computational effort needed in typical S4MPLE docking simulations

Results and discussion

The structures of 3000 selected compounds have been transmitted to our Enamine collaborators. They were able to test 2992 compounds, and the experiments confirmed 29 hits. While this result is objectively low, it is still 2.6 times better than the hit rate found in the random screening of 3200 compounds under identical conditions [START_REF] Borysko | Straightforward hit identification approach in fragment-based discovery of bromodomain-containing protein 4 (BRD4) inhibitors[END_REF]. However, it has to be mentioned that the applied classification models have been trained using publically available SAR data on IC50 values. In contrast, our collaborators have measured the Thermal Shift Assay using Differential Scanning Fluorimetry (DSF). DSF is a biophysical method based on detecting the shift in protein denaturation temperature upon ligand binding, as reported by fluorescent dye interacting with the protein core exposed by heat denaturation. The method is a simple, label-free HTS technology applicable to most soluble proteins, irrespectively of their functions and activities. An in-depth analysis was performed in order to understand:

• How the public-data affinity values used for model building relate to the experimental hit detection criterion (∆Tm) used in DSF?

The original dose-response (such as IC50) activity scores from public sources were shown to be per se rather poorly correlated to the hit selection criterion DSF-ΔTm. The fact that they were not used as such for model training, but first underwent conversion into a categorical variable has most likely had a negative impact on model performance.

• Which of the used models are better at selecting the 29 confirmed hits?

Seventeen hits have been ranked #1 by at least one of the GTM models, while two were ranked #1 by SVM. The other ten hits were selected because of "broader" consensual selection by multiple models that ranked them within the top of the list.

More detailed description of the obtained results is given in our article published in Eur.J.Med.Chem., see below.

Conclusions

A collection of 2M compounds have been subjected to a virtual screening funnel involving classification SVM and GTM models as well as ligand-based pharmacophores trained on publicly-accessible SAR data on BRD4 IC50/pKi from Reaxis and ChEMBL.

Each model has provided a ranked list of 2M candidates according to their likelihood to be active. The consensus application of these approaches has been used to obtain a subset of 12k compounds that have been submitted to a docking procedure. The docking has been used for a further selection of the "best" 3k compounds out of the previously selected pool of candidates. These 3k compounds have been experimentally screened by the Enamine partner using the Thermal Shift Assay method.

Twenty-nine confirmed hits had been detected, which represents 1% of the 3k selected candidates. While the obtained hit rate is still 2.6 times better than the hit rate found in random screening under identical conditions, it is still objectively low. An in-depth analysis of the quality of the used data for the models' training has been performed, as well as the correlation between IC50/pKi and DSF-∆Tm. First of all, it has been shown that public data from different sources cannot be fused into a single and rigorously defined dataset adapted for QSAR modeling. Moreover, it has been shown that the dose-response activity values reported in publicly available databases are weekly correlated with DSF-∆Tm. Last but not least, the retrospective hit analysis has shown that GTM models have outperformed SVM and ligand-based pharmacophores in terms of hits identification.

Combinatorial Libraries of imines

The reversible combination of molecular building blocks via covalent or non-covalent bonds is a cornerstone of the Dynamic Combinatorial Chemistry [START_REF] Lehn | Dynamic combinatorial chemistry and virtual combinatorial libraries[END_REF]. The reversible nature of the reactions between the building blocks leads to the fact that their thermodynamic stability dictates product distribution in the mixture. Once a Dynamic Combinatorial Library (DCL) is exposed to an external effector (like a biological target), it might happen that the latter selectively binds to one (or several) members of the DCL; thus, the equilibrium is shifted according to Le Chatelier principle, which leads to global change of the solution composition. The nature of the effector is not limited to the biological target. It could be physical (like a change of temperature [START_REF] Giuseppone | Protonic and temperature modulation of constituent expression by component selection in a dynamic combinatorial library of imines[END_REF]), or chemical (introduction of an "alien" species to the solution, like a metal ion [85], a protein/enzyme [8]; change of solvent/pH [START_REF] Osypenko | Pattern generation and information transfer through a liquid/liquid interface in 3D constitutional dynamic networks of imine ligands in response to metal cation effectors[END_REF]).

For example, in a hypothetical library, only containing two pairs of building blocks, the outcome is easily predictable -the solution composition will reflect the relative stability of the reactions' products. Once an effector is added, the library composition changes as a function of products affinity to the effector. In the simplest case, only one reaction product selective binds the given effector, leading to the shift of all equilibria in DCL in favor of the formation of the complexed species (Figure 6.1).

each other, more A1B1 is formed, less non-reacted A1 and B1 building blocks remain in solution, therefore the concentration of A2B2 increases. When it comes to the products that are in an antagonistic relationship, it is the opposite. More A1B1 is formed, less A1B2 and A2B1 are present in the solution because of the lack of needed building blocks.

Although the number of known reversible reactions is relatively high, not all of them can be carried out in aqueous media, which prevents their usage in protein-templated DCL.

Moreover, the functional groups of building blocks are not supposed to react with the target themself. The non-exhaustive list of reversible and biocompatible reactions includes imine, hydrazine and acylhydrazone formation; alkene cross-metathesis; disulfide, thioether and boronate ester formation. Imine formation reaction was the first reaction applied to a DCL in the presence of bovine carbonic anhydrase II [8] as a receptor.

Usually, when one wants to use a DCL for the identification of a "best binder," the reaction products are expected to be in almost equal concentrations, since in the case of the biased library, where one or a few constituents would be highly favored, the preferred interaction of a minor constituent with the target may not be strong enough to overturn the equilibrium situation. In some cases, an experienced chemist could guesstimate what reactants should be chosen for the given DCL; however, this approach might be inefficient.

The usage of the software which can estimate speciation, i.e., equilibrium concentrations of all species in solution. This requires a knowledge of equilibrium constants, which can be problematic because of the thermodynamic data availability.

Chemoinformatics models predicting equilibrium constants could be a reasonable solution to assess speciation for any DCL with or without effector. In this work, DCL based on the reaction of imine formation is modeled with and without an external chemical effector. As a tribute to the seminal work [START_REF] Bishop | GTM: The generative topographic mapping[END_REF], human carbonic anhydrase II (CA II) was chosen as an effector to model the adaptive behavior of the imine-containing DCL. The project workflow (Figure 6.4) involves several steps. At the first step, a predictive model for the logarithm of imine formation constant (logK) as a function of the structure is built using experimental data measured in chloroform solution. In the second step, a model for the logarithm of the binding constant (pKi) of organic molecules to human CA II should be prepared on experimental data extracted from the ChEMBL database. Since the latter were measured in water, one needs to scale the imine formation constant determined in chloroform to those in water solution. Once both types of models are available, predicted stability and binding constants obtained could be used as input to a speciation software.

One of the important tasks was to select a representative training set for the imine formation constants modeling. This work is described in section 6.1, followed by the modeling part described in section 6.2. Reported case studies concerned empirically designed DCLs involving a relatively small number of reactants. Our goal is to develop a theoretical (in silico) approach allowing one to predict the species concentration for DCL of any size in the presence or in the absence of effector. The workflow of such in silico speciation of DCL is given on Figure 6.4. It involves two essential steps: (i) selection of a diverse library of imines which serves as a training set in model building, and (ii) preparation of statistical models able to predict equilibrium constants of imines formation and binding constants of protein-ligand complexes in the solution used, in turn, as an input in a speciation software. These two steps of the project are described in two separate sections below.

Application of GTM for diverse library selection

Introduction

The identification of representative and diverse subsets in large libraries of compounds is crucial to medicinal chemistry since a diverse subset of compounds provides more chances to contain a compound with the required type of activity during screening tests. When one or several compounds from a diverse subset have been proven to show a certain level of activity, then a focused library of compounds is being selected. Focused library design implies the selection of similar compounds to the known "hits". Traditionally diverse subsets have been created by having a medicinal chemist select compounds manually based on a series of 2D structures [START_REF] Clark | OptiSim: an extended dissimilarity selection method for finding diverse representative subsets[END_REF]. Although this approach could be successfully used on relatively small datasets, it would become extremely difficult (or probably impossible) when the datasets have thousands and hundreds of thousands of compounds. Moreover, the level of "representativity and diversity" of the manually selected subset could vary from chemist to chemist. The selection of a diverse library from a large dataset of compounds can be approached in several ways: using clustering, using dissimilarity-based methods or using cell-based methods. Clustering implies that the compounds are grouped according to some similarity measure. Then from each cluster, a random compound is selected, thus leading to a library containing a representative of each cluster. Dissimilarity-based methods are relying on the calculation of pairwise distances for every compound of the initial dataset, which is followed by a one-by-one selection of compounds in the diverse library according to a pre-defined rule. The main feature of cellbased selection methods is that they do not require a pairwise calculation of the distances for all the compounds. To provide an efficient application of a cell-based approach, the dimensionality of chemical space defined by the molecular descriptors is usually reduced to a 2D map; the map is divided into zones (cells), and from each cell, a compound is extracted.

Although the diverse libraries are usually used in medicinal chemistry, one can use a diverse subset of compounds as a training set in structure-activity modeling. Diversity is not a fundamental, objective property of a compound collection, but it is a rather practical, problem-dependent and therefore vaguely defined concept. Usually, compound diversity is directly related to their dissimilarity, hence to quantitatively measure the diversity of selected compounds is the most straightforward to encode into chemoinformatics software.

If molecules are objects in the descriptor space, then dissimilarity is directly related to the distances separating them. Depending on the descriptor space and the therein employed metrics, the distances may vary. At first sight, the above seems like a rigorous mathematical basis for compound selection: the degree of dissimilarity is typically illustrated by the distribution histogram of pairwise distances between the compounds.

Despite the large number of different applications of GTM, it has never been used for diverse compounds library design. Since GTM is producing a 2D map, this map could be used in a cell-based approach for diverse library selection. In this context, GTM could be compared to its non-probabilistic predecessor -Kohonen Self Organizing Map (SOM). It has been reported that SOM is a useful tool for focused library design and combinatorial libraries [START_REF] Schneider | Ligand-Based Combinatorial Design of Selective Purinergic Receptor (A2A) Antagonists Using Self-Organizing Maps[END_REF][START_REF] Schneider | Self-organizing maps in drug discovery: compound library design, scaffold-hopping, repurposing[END_REF]. Nettekoven and Schneider [START_REF] Schneider | Ligand-Based Combinatorial Design of Selective Purinergic Receptor (A2A) Antagonists Using Self-Organizing Maps[END_REF] used SOM for the focused-library combinatorial design of selective purinergic receptor (A2A) antagonists. Selzer and Ertl [START_REF] Selzer | Applications of self-organizing neural networks in virtual screening and diversity selection[END_REF] used SOM to select a representative subset of 5000 compounds from a collection of combinatorial libraries containing nearly 100000 compounds in total. Unfortunately, in that work, the method performance in diverse library selection has been reported in a qualitative way as "low/medium/high" diversity of the selected library, which does not allow quantitative comparison of the results obtained with GTM to the results obtained by SOM in a similar task. In this section, the performance of GTM in diverse library selection has been compared to a classical algorithm of diverse library selection -MaxMin [START_REF] Holliday | Definitions of" dissimilarity" for dissimilaritybased compound selection[END_REF][START_REF] Holliday | A fast algorithm for selecting sets of dissimilar molecules from large chemical databases[END_REF].

Data and methods

A substructural search of primary amines and aromatic aldehydes has been done on SciFinder; the results have been sorted by the number of citations, putting the most cited compounds on the top of the list. The reactants have been selected according to the following criteria:

• Every compound should contain one single amine/aldehyde group, thus leading to only one possible reaction product.

• In the case of amino-acids, the COOH group has been changed into COOMe.

• Thiol containing compounds have been rejected.

• Reactants having a molecular weight > 400 g/mol have not been taken into account.

• Long-chained acetals -C(OR1)(OR2) have been changed to -C(OMe)(OMe).

Two datasets of reactants have been obtained: the primary amines dataset containing 300 molecules and 400 compounds dataset of aromatic aldehydes. Interaction between all the selected aldehydes and amines could give 120 000 imines. Since MaxMin is a very timeconsuming method of O(n 2 N) complexity (n and N are the numbers of objects in the initial set and diverse subset, respectively), we decided to perform some methodological tests on a smaller set containing 42658 imines resulted from reactions of 154 amines and 277 aldehydes. Five different ISIDA fragmentations have been used (see section 6.1.2.2 for details):

Fragmentation scheme Meaning

IA-FF-FC-2-3 Sequences of atoms colored by force field properties and formal charge with a length from 2 to 3 IAB-FC-1-6 Sequences of atoms and bonds with a formal charge on atoms, with a length from 1 to 6 IIRAB-FF-1-3 Circular fragments of atoms and bonds with restricted length from 1 to 3 atoms, colored by the force field

IAB-FC-2-4
Sequences of atoms and bonds with a formal charge on atoms, with a length from 2 to 4

IIRA-FF-FC-1-2
Circular fragments of atoms with restricted length from 1 to 2 colored by the force field and containing information about the formal charge

Dissimilarity-based methods

Since the task is to select the most dissimilar compounds, one should select the compounds that are most different from the ones already selected in the formed set. As the metric of dissimilarity, the Soergel distance [START_REF] Gower | Measures of similarity, dissimilarity and distance[END_REF][START_REF] Willett | Chemical similarity searching[END_REF] has been taken. The following formula defines Soergel distance:

= 1 - ∑ ∑ + ∑ -∑
Where xjA and xjB represent descriptors vectors of compounds A and B, respectively. It has been shown that Soergel distance can be used as a metric only when the values of the descriptors are non-negative [START_REF] Lipkus | A proof of the triangle inequality for the Tanimoto distance[END_REF]; therefore the ISIDA descriptors are well adapted for this study.

Holliday et al. [START_REF] Holliday | Definitions of" dissimilarity" for dissimilaritybased compound selection[END_REF] evaluated the performance of 4 dissimilarity-based methods -MaxMin, MaxMax, MaxSum and MaxMed. These methods follow similar library selection algorithms:

1. Calculation of all N(N-1)/2 pairwise distances (dissimilarities) for N compounds in the dataset.

2. Selection of random compound (Compound 1) from the initial library and addition of it to the subset S.

3. Compound 2 is the most remote compound with respect to Compound 1.

4. Identification of the most dissimilar compound to the already selected compounds from the initial dataset using a distance-related score and its addition to the subset S.

Repetition n-2 times of step 4.

The difference between MaxMin, MaxMax, MaxSum and MaxMed concerns the rules used to evaluate a score used to identify the next object to be included in the subset S.

Let Dij be the dissimilarity between the i-th molecule in the initial dataset and j-th molecule in the subset S. As a function of the algorithm, and the following scores are used: MIN{Dij} for MaxMin, MAX{Dij} for MaxMax, Σ{Dij} for MaxSum and Median{Dij} for MaxMed.

The object having a maximal score is added to the subset S.

It has been shown [START_REF] Holliday | Definitions of" dissimilarity" for dissimilaritybased compound selection[END_REF] that MaxMax and MaxMed often led to subsets containing too similar compounds. MaxSum method preferentially selected compounds that were located "at the corners" of the chemical space. These drawbacks are minimized within the MaxMin method, which ensures a relatively uniform selection of compounds from different areas of the chemical space. The most significant drawback of all dissimilarity-based methods is the necessity to calculate the distance matrix for all the compounds. The complexity is proportional to O(n 2 N) [START_REF] Holliday | A fast algorithm for selecting sets of dissimilar molecules from large chemical databases[END_REF], where N is the number of compounds already selected, and n is the total number of compounds in the initial library.

Here, the MaxMin algorithm was applied to both reactants and products. For a product-based approach, the MaxMin algorithm has been directly applied, and 225 imines have been selected. Since the algorithm takes the first compound randomly, 100 libraries have been selected, each having a different seed. In the case of the reactant-based approach 15 aldehydes and 15 amines were selected, and the corresponding imines to the pairwise reactions between the selected reactants were directly extracted from the imine library.

Since the number of aldehydes and amines is rather low (277 aldehydes and 154 amines)

and the algorithm is dependent on the choice of the first compound in the list, for the selection of the most diverse library of amines and aldehydes, an "exhaustive MaxMin" approach has been applied, where each compound has been used as the seed. Thus, the usage of the "exhaustive MaxMin" for the reactant-based approach ensures that the selected libraries of reactants are indeed the most diverse. To compare MaxMin applied on products and on reactants, ten the most diverse libraries of aldehydes and ten most diverse libraries of amines have been selected to generate 100 most diverse libraries of imines.

Cell-based diverse-library selection using GTM

In order to use GTM as a method for diverse library design, the compounds were projected on the map, and the map itself was virtually evenly divided into n cells, where n equals to the number of needed compounds (therefore in this study n=225). From each cell a random compound has been extracted. GTM parameters have to be optimized for this task in order to obtain a map covered as evenly as possible by all the compounds of the initial dataset. To do so, the same protocol of parameter optimization involving the Genetic Algorithm (GA) [START_REF] Horvath | An evolutionary optimizer of libsvm models[END_REF] described previously has been applied. In these calculations, normalized Shannon entropy has been used as the scoring function. The following formula calculates Shannon entropy [START_REF] Lin | Parallel Generative Topographic Mapping: an Efficient Approach for Big Data Handling[END_REF]95]:

= - CumR log(CumR )
Where CumRk is the cumulated responsibilities of compounds in the node k. However, in the GA optimization, the normalized Shannon entropy has been used:

= log ( ) × 100
Where N is the total number of nodes. The normalized entropy ranges within [0; 100],

where 0 means that all the compounds are mapped into the same node, and 100 means that the compounds are covering the map uniformly. In other words, the Shannon entropy corresponds to the homogeneity of the distribution of the source library over the map area.

The higher the entropy is, the more the objects are dispersed on the map. In such a way, five "best" different descriptors spaces leading to "optimal" maps have been selected

(Table 6-1).

Performance evaluation

The performance of the cell-based approach using GTM has been compared to the performance of the MaxMin algorithm, and the library of randomly selected compounds has been used as the baseline. The quality of selected libraries has been considered taking into account the following criteria:

• Diversity of the selected library. Here two diversity scores were used "All-Soergel" (AllS) and "Min-Soergel" (MinS); the first calculates the average of pairwise distances for all selected compounds in subset S, whereas the latter calculates the average of the distances to the closest neighbor of each compound in the selected subset S.

= 1 -1 = MIN , ,
Notice that the "All-Soergel" score accounts for all pairwise distances in the selected library, whereas "Min-Soergel" shows how dissimilar are the closest neighbors.

• Data coverage. A diverse library of 225 imines will be used to build a model for logK. Which, in turn, will be applied to assess this thermodynamic parameter for the initial set of 42658 imines. However, because of the fragment control applicability domain (see section 3.1), some predictions are considered unreliable. In such a way, data coverage is defined as a ration of the molecules for which predictions are considered reliable to the size of the initial data set.

Results

According to the Min-Soergel score, MaxMin applied to the dataset of products selects more diverse libraries than any other studied approaches (Figure 6.6, Table 6-2). In contrast, MaxMin applied to reactants is much less diverse because any imine in the selected library shares a common aldehyde-substructure with 14 imines and a common amine-substructure with another 14 imines, which reduces the distance to the closest neighbor. According to both All-Soergel and Min-Soergel scores, libraries obtained by cellbased selection are less diverse compared to those selected with the MaxMin algorithm. It can be explained by the fact that MaxMin maximizes the diversity of the selected subset explicitly, while cell-based methods do not consider the distances between the compounds in chemical space. Moreover, MaxMin mainly selects the objects on the "border" of chemical space, while GTM focuses on dense clusters of compounds, and hence, the objects remote from the manifold are not be adequately taken into consideration. On the other hand, the diverse libraries obtained by cell-based selection have several advantages over the classical dissimilarity-based diverse library selection method. First of 123 all, this concerns high data coverage (Figure 6.7, Table 6-3). Since every zone of the 2D latent space corresponds to the related zone in the initial N-dimensional chemical space, regular selection objects from the "cells" on GTM corresponds to regular sampling from the initial space. The second advantage of the cell-based approach over MaxMin is the high speed of calculations (Figure 6.8). The slowest step of the MaxMin algorithm is pairwise distance calculation. For instance, this step took 1h on average for the calculation of pairwise distances of 42658 compounds. On the other hand, in a cell-based approach, even selection objects on a 2-dimensional map takes several seconds. Moreover, in the latter case, the time of calculations does not depend on the size of the initial dataset One can also envisage a two-step workflow combining MaxMin and cell-based algorithms. In the first step, the cell-based method is applied in order to select an "intermediate" diverse library whose size is larger than of the final set (225 compounds).

This "intermediate" set serves as a source library for the MaxMin algorithm, which will select 225 the most dissimilar compounds. The results show (Figure 6.9) that the combination of two approaches is a reasonable trade-off between the diversity of the selected library and the time of computations.

Discussion

The above results demonstrate that GTM is an acceptable method to select diverse library: it performs reasonably well, it is very fast, and its data coverage is almost 100%.

On the other hand, all 225 imines selected with GTM may contain unique reactants. This may be a problem for the budget of the experimental laboratory, which needs to purchase 225 amines and 225 aldehydes.

In this regard, the dissimilarity-based MaxMin algorithm applied to reactants might be a reasonable solution despite its relatively small data coverage. Therefore, we decided to use the latter for the selection of a training set for logK of imines formation modeling in a two-steps procedure. In the first step, five different diverse libraries larger than 225 compounds were selected using descriptors spaces mentioned in Table 6-1. Their overlap resulted in 15 aldehydes and 15 amines present in all five individual libraries. Then, because of the recommendation of experimentalists, the subset of aldehydes was extended to 24 aldehydes. Their pairwise combinations result in a library of 360 imines, out of which 276 imines were synthesized, and their logK were measured using NMR spectroscopy (see experimental details section).

Chemoinformatics driven assessment of speciation in dynamic combinatorial libraries

Structures related to logK and pKi datasets were standardized following the procedure implemented on the virtual screening server of the Laboratory of Chemoinformatics at the University of Strasbourg (infochimie.u-strasbg.fr/webserv/VSEngine.html) using the ChemAxon Standardizer. The evolutionary model tuning of the Support Vector Regression [START_REF] Horvath | An evolutionary optimizer of libsvm models[END_REF] approach was applied to grow both stability and CA affinity Support Vector Regression (SVR) models. As the approach can select the best suited molecular descriptors out of a user-provided pool of potentially useful descriptor sets, in both cases, the optimizer was given the freedom to choose its preferred descriptor spaces. However, given the different scopes of the approaches (stability -based on imines representing a combinatorial core of the envisaged DCL/affinity -based on public compounds with reported CA affinity data, most of them not being imines), distinct pools of ISIDA fragment descriptors were fed as possible input. For stability, a series of 18 customized ISIDA fragmentation schemes have been selected as a result of a genetic algorithm; then, in each of the 18 descriptor spaces, an additional evolutionary model tuning has been done. The protein binding predictor was grown, starting with the "default" pool of 100 ISIDA descriptors schemes Note that each data point is labeled either "exact" or "estimated". These labels come from NMR limitations. The "estimated" label was given to those logK of imines when the concentration of products/reactants were hard to identify precisely.

Experimental details

NMR measurements. Imines were prepared directly in 5 mm NMR tubes by mixing 200 mM stock solutions of components (60 µL each) and diluted with 480 µL of CDCl3 to reach the final concentration of imines of 20 mM. The NMR measurements were performed after 24 h of equilibration at room temperature. 1 H MNR spectra were recorded on 500 MHz Bruker spectrometer with an automated sampler, using standard parameters.

Solvent. In this study, the deuterated chloroform CDCl3 was used to perform all the tests of imine formation for several reasons. Namely, (i) many organic molecules are soluble in chloroform; (ii) the formation of imines in this solvent, in general, is quantitative;

(iii) it is a solvent of choice for routine NMR. Before use, chloroform was filtered through the pad of basic alumina to remove residual acid always present. Next, several milliliters of Mili-Q water were added to the bottle to obtain a saturated solution of water in chloroform. This strategy helps to ensure control of the water content to be relatively constant. During the imine formation, one molecule of water is produced, and it participates in the equilibrium. However, it is pretty difficult to measure precisely the amount of water by NMR, so saturation of chloroform with water should solve this problem.

The concentrations of imines, as well as of reactants, were, thus, obtained by integrating the corresponding signals. For most of the cases, the concentrations have been easily identified; hence the associated logK was precisely calculated. However, in some "extreme" cases, the signals' intensity was so weak (either the concentrations of both reactants or imine being low), that it was nearly impossible to integrate the peaks correctly.

In this case, the obtained value of logK was labeled as "estimated".

logK modeling in chloroform

16 SVR individual models, each built in 16 selected descriptor spaces (Table 6- they, in 90% they produce imines with a negative logK. As a counterpart to "inert," some "reactive" compounds have been found, in 80% of the cases when a "reactive" aldehyde or amine is involved in imine formation reaction, the outcoming logK values are usually high (> 4). It should be noted that "inert"-amines often bear methyl ester fragment, while 3,5dibromobenzaldehydes are "inert" (Figure 6.15). On the other hand, "reactive" reactants share a common substructure of 3,4,5-trimetoxibenzene and benzodioxole (Figure 6.16).

Last but not least, some "versatile" reactants were also identified. These reactants can be potentially involved in any imine formation since related logK values range from -8 to 8.

The structures of "versatile" reactants are shown in the Figure 6.17. constant of 7.5). The application of the equations mentioned above resulted in ∆K chl->w ≈ 0.5 log units. Note that in water, the solvent is a reaction product. Hence, the massive presence of water would displace equilibrium towards hydrolysis, and might also "capture" a part of the aldehyde reagent under the form of hydrates. However, at this stage, there is no experimental data available to support our calculations, which limits us to the debatable but not absurd working hypothesis that the relative order of imine stability is not too much affected.

Models description

The affinity model for human CA II is publicly available within the property prediction tool on the web server of the laboratory of Chemoinformatics (see section 6.3). It is composed of a consensus predictor based on the top 5 evolved models using the top 5 best suited ISIDA descriptor spaces. Each of the individual models (Table 6-6) is applied to the compounds submitted to the webserver. However, the output of the ones containing the compound to predict within its specific AD (according to the descriptor-specific fragment Control rules) is preferentially used to calculate the returned consensus (mean) value. The mean of all predictions, irrespectively of AD compliance, is also returned -it may be used as a low-trust estimator for compounds that are out of the AD of all the five individual models. [a] Fragmentation scheme nomenclature in column 1 denotes the fragment type (I-sequence, II-circular fragments), the nature of captured information (A-atom types are captured, B -bond orders are captured), the coloring scheme (FF -force field type-based labeling supersedes default labeling by atomic symbol), other options (P -atom pair counts only, FC -formal charges are considered).

The developed pKi model has been applied to a set of 80400 imines within the applicability domain of the SVR consensus model for imines equilibrium constant. The distribution of predicted pKi values shows that most of the imines should have a pKi value between 5 and 6 (see Figure 6.21).

Speciation assessment

The ChemEqui program [START_REF] Solov'ev | Supramolecular complexes: Determination of stability constants on the basis of various experimental methods[END_REF] has been used in this project. For the prediction of equilibrium concentration, it uses two groups of equations:

• The law of mass action (Brinkley's representation), where Ci is a reaction product i, βi is the formation constant of component i νij is the stoichiometric coefficient in the ith chemical equilibrium involving basic component j with concentration Cj and m is the number of basic components.

= +

• The law of mass conservation, where C 0 j and C Σ j are respectively the initial concentration and analytical concentration of component j, r is the number of reactants.

= =

Then, the following function is minimized: only logK varied in the range of 2.5 -4.25. The results (Table 6789) have shown that the selective binding with the selected imine takes place if:

• logK of selected imine drops by 1.25 with respect to its initial value.

• logK of selected imine and its "agonistic pair" changes by 0.6 • logK of selected imine decreases/increases and its "antagonistic pair" increases/decreases by 0.75 It follows from the above computational tests that

• logK of considered imines should be as close as in order to obtain their quasiequivalent distribution.

• pKi values should differ by at least 1.5 log units.

Following these recommendations, four imines resulted from the interaction between two aldehydes, and two amines shown in Figure 6.19 have been selected. Their predicted logKwat and pKi values are given in Table 6-10.

Models implementation 6.3.1 Predictive models of logK of imine formation in chloroform

The obtained consensus SVR models were uploaded on the Predictor service of the laboratory of Chemoinformatics (http://infochim.u-strasbg.fr/cgi-bin/predictor_dcl.cgi). A short "user-guide" (Google Chrome v. 84.0 or Mozilla Firefox v. 78.0 browsers) for proper usage of the implemented models is given below:

1. First, select the "DCL model" on the left menu of the predictor and set the "Select a general kind of property" to DCL and "Select a property to model" to Imine_Formation_LogKeq_Reg.

2. The user can either draw the imine for which he wants to receive a prediction or, in the case of multiple imines, he can give an sdf. Note that the Predictor tool can treat only 100 compounds. Once one of the two input possibilities were chosen, click on the "Submit" button.

Predictive models of pKi of human CA II

These models have been uploaded on the web server of the laboratory of Chemoinformatics (http://infochim.u-strasbg.fr/webserv/VSEngine.html) with the help of Dr. Dragos Horvath. A short "user-guide" (Google Chrome v. 84.0 or Mozilla Firefox v.

78.0 browsers) for proper usage of the implemented models is given below:

1. Enter a preferred username and password in the top left corner of the webpage in order to connect to the web server.

2. Once connected to the web server, click on "QSAR-based Property

Predictions" located on the top left corner of the webpage.

3. The webpage will be refreshed and a new field will appear. asking the user to enter the project name.

4. On a new webpage, the user can either draw a compound or provide an sdf or smiles file containing multiple compounds. The web server will then standardize the provided compound(s). Do note that for large files having thousands of compounds, it may take some time.

5. Once the standardization is done, the user will be redirected to a webpage of all the implemented models on it. The model concerning the current project is called CarbAnhydrII-pKi-CHEMBL205. Once the model is selected, the user can proceed to the predictions by clicking on the "GO!" button. The webserver tool will start the needed descriptors generation and eventually will proceed to the predictions.

6. Once the predictions are made, the webserver tool will generate a .csv file containing the predicted values.

7. The obtained file will contain a short description of the used model as well as the number of the compounds that have been sent to the prediction tool. An explanation for every column will also be given.

Conclusions

To understand the behavior of a DCL, the relative propensity of the formation of DCL products and their relative affinities for the target are needed. The first steps in this direction have been done, including the preparation of predictive models for imines stability constants and affinities of organic molecules for the human CA II. While this work underlies the conceptual workflow and reports useful and publicly available models, experimental proof of the reported here in silico approach is needed.

The calculation of speciation of any solution does not present any technical difficulty due to numerous software; its utilization remains constrained by the availability of experimental data. The study presented herein aimed to overcome that constrain by the usage of current chemoinformatics methods and tools. First, from the pool of 400 most cited aromatic aldehydes and 300 most cited primary amines (according to SciFinder), diverse subsets of 24 aldehydes and 15 amines (resulting in 276 imines) have been selected.

These imines have been synthesized. Their logK measured and the obtained data served as a training set for building SVM classification and regression models. These models showed high "extrapolation potential" since they were able to provide reliable predictions for more than 80000 imines from the initial pool of 120000.

Although the models themselves are useful and the obtained predictions can be already used for the speciation of a DCL without any effector, they are not sufficient enough for the modeling of a DCL in the presence of any effector. In order to model the presence of the effector, it is needed to quantitatively know the affinity of each constituent in solution to it for this purpose. ChEMBL database has been used as a source of binding affinity (pKi) on human Carbonic Anhydrase II. More than 4300 compounds have been extracted and used to train SVM regression models.

Developed in this work, predictive models for logK and pKi allowed us to overcome the "experimental data constrain" and thus make possible the usage of speciation software.

In considered here hypothetical DCL, the selected quartet of imines has close predicted values of logK, thus ensuring a quasi-equal concentration of imines. Moreover, for imines, related pKi values differ by at least 1.5 log units, thus ensuring a certain level of selectivity.

It follows that on the example of "imine-based" DCL modeling, the available speciation software coupled to chemoinformatics tools could be, in principle, used to any type of DCL and any chemical/biological effector if there is a sufficient amount of data for the modeling.

Conclusion and Perspectives

This work was devoted in the first place to the application and study of Generative Topographic Mapping in virtual screening. GTM has also been used in a task of selection of diverse libraries of imines, and the quality of the obtained libraries has been compared to the classic dissimilarity-based method of diverse library selection -MaxMin. Last but not least, the first steps towards the modeling of protein-templated dynamic combinatorial libraries have been made.

Two projects focused on the application and exploration of GTM have been done.

The first project was rather methodologically-oriented since it provided several universal maps (multi-target GTM-based classification models) that can accommodate more than 1.5M compounds extracted from ChEMBL and to discriminate actives from inactives with high accuracy for 617 targets out of 618. Moreover, this study has shown that the maps are complimentary since each map has been built in different descriptor spaces; therefore, if certain target-specific activities are poorly predicted by one map, there would be another map that will be able to do it. The usage of DUD targets as external-validation sets helped to identify that the correlation between predictions quality in cross-validation is weak. This fact motivated the usage of the universal maps in the consensus model since one could not tell apriori what map will show better results in "real-life tasks". It has been shown that universal maps applied in consensus provide undeniable advantages such as i) 100% of data coverage for most of the targets; ii) higher performance in cross-validation and external validation according to BA and ROC AUC; iii) higher enrichment factor for top 100

predicted "active" compounds.

The goal of another project was to find new inhibitors of Bromodomain 4 by virtually screening a collection of 2M compounds. In this project, a virtual screening funnel was composed, including the building of ligand-based pharmacophore models, SVM and GTMbased classification models on publicly available data extracted from REAXIS and ChEMBL. The obtained models have been used in consensus and selected 12k compounds that have been predicted by most of the models to have the highest probability of being active. Then, this subset of 12k compounds has been subjected to a docking procedure, which selected 3k most potent compounds. These compounds have been tested by our collaborators from Enamine. Out of 2992 tested compounds, only 29 compounds have been identified as "active". While this result is still 2.6 times better than the screening of a set of 3k randomly selected compounds, 29 confirmed hits objectively is a low success rate. First of all, it has been shown that public data from different sources cannot be fused into a single and rigorously defined dataset adapted for QSAR modeling. It has been found that the active/inactive labels of training data have been assigned according to pKi and IC50 values, while the experiment has been done using Differential Scanning Fluorimetry -a method that measures ∆Tm (thermal denaturation temperature). Moreover, it has been shown that the correlation between IC50/pKi and ∆Tm. is very low. Still, it has been shown that GTM classification models were able to find 24 out of 29 confirmed hits.

For the first time, GTM has been used for a selection of a diverse library of compounds. The quality of the selected libraries has been compared to the quality of libraries obtained by MaxMin -a classic dissimilarity-based method for diverse libraries selection. Since the term and the metrics of "diversity" are vaguely defined, a score based on Soergel distances of the compounds included in the diverse library has been used. It has been shown that GTM as an individual model cannot provide diverse libraries with the same level of dissimilarity as MaxMin; however, GTM-based diverse library selection selects a diverse library that is more "representative" than the libraries selected by MaxMin.

Moreover, it was found that the application of GTM for a pre-selection of a bigger diverse library ("intermediate" diverse library) followed by the application of the MaxMin method on the "intermediate" diverse library gives much better results than the application of GTM individually.

The first steps were made towards the modeling of the dynamic combinatorial libraries. To understand the behavior of a DCL, the relative propensity of the formation of DCL products and their relative affinities for the target are needed. The calculation of speciation of any solution does not present any technical difficulty due to numerous software, but its utilization remains constrained by the availability of experimental data.

First, from the pool of 400 most cited aromatic aldehydes and 300 most cited primary amines (according to SciFinder), diverse subsets of 24 aldehydes and 15 amines (resulting in 276 imines) have been selected. These imines have been synthesized, their logK measured, and the obtained data served as a training set for building SVM classification and regression models. The obtained models have shown high data coverage of the initial pool of 120k imines (67% of the data are in the AD of the consensus model) as well as high predictive performance. In order to model the presence of the effector, it is needed to quantitatively know the affinity of each constituent in solution to it. For this purpose, the ChEMBL database has been used as a source of binding affinity (pKi) on human Carbonic Anhydrase II. More than 4300 compounds have been extracted and used to train SVM regression models. Both regression models allow the overcoming of the initial "experimental data constrain" and thus make possible the usage of speciation software. The mutual usage of logK and pKi regression models, as well as ChemEqui speciation software, lead to a possibility of modeling a hypothetical DCL containing 2 aldehydes, 2 amines and human CA II as an effector.

Perspectives

The candidates for DCL in the presence of biological target should be selected in such a way that (i) neither aldehydes nor amines efficiently interact with the effector and, (ii) only one imine firmly binds to the protein. For this purpose, a series of ligand-to-protein docking calculations could be envisaged to gain a microscopic insight into imine-protein interactions.

Experimental measurements of DCL equilibria in aqueos solution are very welcome.

They would help to build new predictive model for logK in water and to prove our suggested here protocol of logK rescaling from one solvent to another one.

Iuri CASCIUC

Cartographie chimique et modélisation de systèmes complexes

Résumé

Cette thèse concerne l'utilisation de la Cartographie Topographique Générative (Generative Topographic Mapping -GTM) pour l'analyse et la visualisation de jeux de données, le criblage virtuel et la conception de chimiothèques. La performance en criblage virtuel de modèles GTM de classification multi-cibles (uGTM) a été étudiée et l'utilisation de plusieurs uGTM en consensus a été proposée. Un criblage virtuel utilisant une combinaison de la GTM avec d'autres techniques de chémoinformatique a permis de découvrir 29 nouveaux inhibiteurs de BRD4 dont l'activité a été prouvée expérimentalement. La GTM a été comparée à la méthode MaxMin comme outil de conception de chimiothèques. Il a été trouvé que malgré le fait que les chimiothèques obtenues avec MaxMin sont plus diverses que celles obtenues avec la GTM, cette dernière est plus rapide et peut être appliquée à des jeux de données plus larges. Un protocole de modélisation pour l'analyse de spéciation de bibliothèques combinatoires dynamiques basées sur la réaction de formation d'imines en absence et en présence d'une protéine a été proposé. Les modèles développés sont disponibles au public sur le site de Laboratoire de Chémoinformatique.

Mots-clés : GTM, QSAR ; criblage virtuel, visualisation de données, conception de chimiothèques, Bibliothèques Combinatoires Dynamiques.

Résumé en anglais

This work concerns application of Generative Topographic Mapping method to different tasks including data analysis and visualization, virtual screening and library design. Performance of multi-target GTM-based classification models (uGTM) in virtual screening was investigated and consensus usage of several uGTMs has been suggested. Virtual screening involving a combination of GTM with some other chemoinformatics techniques allowed to discover 29 new BRD4 inhibitors, activities of which were experimentally confirmed. As a library design tool, GTM was compared to the MaxMin method. Although diversity of MaxMin libraries is systematically larger than those obtained with GTM, the latter is much faster and, therefore, can be recommended for large datasets. A modeling workflow for speciation analysis in imine-based Dynamic Combinatorial Libraries in absence and presence of a protein has been suggested. Developed models are publicly available at the site of the Laboratory of Chemoinformatics.

Key words: GTM, QSAR, virtual screening, data visualization, library design, Dynamic Combinatorial Libraries

Figure 1 . 1 :

 11 Figure 1.1 : Préparation d'une carte générative topographique (GTM) pour un espace chimique définie par les descripteurs moléculaires

Figure 1 . 2 :

 12 Figure 1.2 : La performance cumulée de cartes exprimée en nombre de cibles ayant un BA supérieur au seuil établi en fonction de nombre des cartes utilisées au criblage

Figure 1 . 3 :

 13 Figure 1.3 : Une des cartes génératives topographiques utilisées lors du VS. Les zones rouges et bleus sont peuplées par des composés actifs et inactifs, respectivement. Les régions ayant des couleurs intermédiaires correspondent aux zones peuplées simultanément par des composés de deux classes.

  proposer une chimiothèque diverse de 225 imines. Deux algorithmes ont été comparés : l'algorithme « traditionnel » MaxMin et un algorithme innovant basé sur la GTM. L'algorithme MaxMin maximise les distances entre les individus de la chimiothèque diverse et sélectionne la molécule la plus éloigné de l'ensemble des molécules déjà choisis. Pour cette raison, son utilisation pour échantillonner une très grande chimiothèque conduit rapidement à des temps de calcul prohibitifs. L'algorithme innovant basé sur une carte GTM exploite les deux dimensions de la carte. En effet, sur un nombre réduit de dimension, une approche directe qui consiste à diviser la carte en cellules de surface égale puis à tirer au hasard des représentants dans chacune d'elles, est très efficace. La notion de distance entre individus dans la sélection est écartée, mais l'échantillon est bien plus représentatif. Ici, la carte était divisée en 225 zones égales, puis dans chaque zone un composé a été extrait aléatoirement. Deux critères de performance ont été utilisés : distance moyenne de Soergel (<S>) et le taux de couverture de données. Le paramètre <S> est calculé comme la distance de Soergel moyenne entre tous les 225 composés sélectionnés. Le taux de couverture désigne le pourcentage de jeu de données initial (42658 imines) qui a un analogue parmi les 225 composés sélectionnés. Une sélection aléatoire de 225 imines a été réalisée pour servir de référence. Les résultats présentés dans la Figure 1.4 montrent que MaxMin permet de sélectionner les bibliothèques plus diverses et, pour peu que le choix ait été réalisé sur les produits de réaction, la sélection offre une excellente couverture du jeu de données, mais cela nécessite beaucoup plus de temps et de ressources informatiques. Toutefois, le taux de couverture des chimiothèques sélectionnées à l'aide de la GTM (98%) est supérieur à celui de MaxMin, mais offre des garanties sur la distance séparant les composés sélectionnés.

Figure 1 . 4 :

 14 Figure 1.4 : La comparaison de la performance de GTM (en vert) en constitution de la bibliothèque diverse des imines avec la méthode traditionnelle (MaxMin) appliqué aux jeux de données de réactant (rouge) et produits (bleu), ainsi qu'avec le choix aléatoire (violet). Le diagramme de gauche représente la distance de Sorgel moyenne, celui de milieu la couverture de l'échantillon et celui de droite la vitesse de constitution d'une bibliothèque sélectionnée.

Figure 1 . 5 :

 15 Figure 1.5 : Spéciation dans un DCL hypothétique de 2 amines et 2 aldéhydes dont la spéciation avant / après l'addition de l'anhydrase carbonique est estimée en utilisant les constantes d'équilibre prédites.

Figure 2 . 1 :

 21 Figure 2.1: Ares of successful applications of GTM.

  molecular weight. An obvious flaw of these descriptors is the impossibility to discriminate between isomers. 2D molecular descriptors are based on a two-dimensional representation of the molecule (a molecular graph). This representation provides the information on atoms connectivity, therefore overcoming the flaw of 1D descriptors. Topological indices and molecular fragments are good examples of this type of descriptors. 3D molecular descriptors are obtained from the 3D structure of the compound. They include the quantitative values obtained by quantum mechanics (such as HOMO/LUMO energies of the compound, dipole moment or electrostatic potential), ovality of the compound and van der Waals volume. In this work, ISIDA [54-56] descriptors are used. They are 2D descriptors that encode a compound structure by counting the number of occurrences of different substructural fragments. These fragments could be linear sequences, augmented atoms

Class 2

 2 False Negative (FN) True Negative (TN)

  behavior. In this thesis, k-fold cross-validation has been used. It consists of the dataset division into k subsets (folds); k-1 folds are used alternatively for model training, and the last fold is used for testing. At the end of this procedure, every molecule has received a prediction exactly once, and these values are used for the calculation of the model's performance parameters (BA for classification, RMSE and/or R 2 for regression).

Figure 3 . 1 :

 31 Figure 3.1: Schematic representation of a 3-fold cross-validation procedure. The initial dataset is divided into three parts; on each fold, the model is trained on two parts, and it is applied to the associated test part. At the end of the procedure, all the predicted values of "test" subsets are used for the model's evaluation.

Figure 3 . 2 :

 32 Figure 3.2: A schematic representation of a separable problem in 2D space. The margin providing the widest separation as well as the hyperplane are defined by the support vectors [65].

  class landscape represents the GTM-based classification model. To obtain a class landscape, a class is attributed to each node of the grid by averaging the responsibilities rkn(Ci) over the number of compounds Nci of the training set that belongs to the i-th class. The conditional probability P(k|Ci) of the new object close a node k is calculated:

Figure 3 . 3 :

 33 Figure 3.3: Examples of three types of landscapes. The GTM has been applied to the dataset containing 6.7k compounds with known activities for vascular endothelial growth factor receptor 2 (CHEMBL279). Class landscape (a) shows the distribution of compounds of two classes: active (red) and inactive (blue). logS (solubility) has been used as a 3 rd axis to build the property landscape (b). Density landscape (c) shows the population of the zones of the map.

Figure 4 . 6 :

 46 Figure 4.6: Activity landscape for ChEMBL301.

. 2 )

 2 of this project involved several machine learning methods (GTM and SVM that have already been described in the Methods chapter), ligand-based pharmacophores and docking. SVM classification models, uGTM BRD4 landscapes, "local" GTM BRD4 landscapes were used in consensus with ligandbased pharmacophore models. Each model has treated the library of 2M compounds, and it ranked them by the likelihood of being active to BRD4, putting the most active compounds on the top of the list. Basing on the predictions of every individual model, 12000 compounds have been selected and sent to the docking procedure. A compound was selected if it was ranked among the top 10 by at least two models, or it was placed among the top 10000 compounds by > 50% of the models.

Figure 5 . 1 :

 51 Figure 5.1: On the left -the structure of Bromodomain -4 α-helices linked by two loops BC and ZA. On the right -BRD inhibitor (Ischemin) interactions with the protein. Binding site residues are shown in sticks. Note that only one hydrogen bond (red dotted line) is made with Aspargine1168 [46].

Figure 5 . 2 :

 52 Figure 5.2: Applied virtual screening protocol.

Figure 5 . 3 :

 53 Figure 5.3: LigandScout pharmacophore features.

Figure 5 . 4 :

 54 Figure 5.4: Confirmed hits projected on one of the used GTM landscapes. The red and blue zones of the map are populated by, respectively, active and inactive compounds. The regions of the map colored in "intermediate" colors are populated by the compounds of both classes.

Figure 6 . 3 :

 63 Figure 6.3: A DCL composed of 2 pairs of building blocks. The pairs in an agonistic relationship are shown with green lines. The compounds being in an antagonistic relationship are shown with red lines.

Figure 6 . 4 :

 64 Figure 6.4: Workflow for in silico DCL speciation.

Figure 6 . 5 :

 65 Figure 6.5: Density landscapes of a set of 42658 imines built in 5 different descriptor spaces that have been used for cell-based diverse library selection. Maps parameters are given, where k is the square root of the number of nodes, and m is the square root of the number of RBFs.

Figure 6 . 6 :

 66 Figure 6.6: Diversity of selected libraries according to the "All-Soergel" score (top) and "Min-Sorgel" score (bottom). Each value is a score average of over 100 diverse libraries selected with different random seeds.

Figure 6 . 8 :

 68 Figure 6.8: The time needed for every approach to select a diverse library of 225 compounds on a computer with Intel(R) Core(TM) i7-6900K CPU and 16 GB of RAM.

Figure 6 . 9 :

 69 Figure 6.9: Diversity of selected library by consecutive application of cell-based approach and MaxMin algorithm according to Min-Soergel score.

Figure 6 . 10 :

 610 Figure 6.10: Model building workflow applied for both logK and pKi.

Figure 6 .

 6 Figure 6.11 shows the distribution of measured logK values for 276 selected imines.

  4),contributed to consensus calculations. The resulting consensus model well performs in 5fold cross-validation: the determination coefficient R²=0.93 and the RMSE=0.62 log units.A plot of experimental vs. predicted logK values is shown in Figure6.12.

6. 2 . 1 . 4

 214 Chemical space analysis of a set of 80400 iminesTHE developed SVR consensus model has been applied to a set of 80400 hypothetical imines identified as being inside of its applicability domain. Distribution of the predicted logK values shows that only 30% of imines have logK> 3 and, therefore, are suited to be used in DCL.GTM has been used to visualize a chemical space of 80400 imines. The manifold has been built on 276 imines with experimentally determined logK values. Then all the imines that are in the applicability domain of the SVR consensus model have been projected on the map, thus creating a property landscape (Figure6.14). One can notice that the majority of the compounds having negative logK values are located on the central and right-central side of the map (green and light/dark blue colors); the compounds having high logK values (orange/red) color are located on the top left and bottom right side of the map. An in-depth analysis of chemotypes has been done in order to see what reactants are usually linked to low and high logK values. Some aldehydes and amines have been identified as "inert" since

Figure 6 . 13 :

 613 Figure 6.13: Distribution of predicted values of logK of imine formation in chloroform.

Figure 6 . 14 :

 614 Figure 6.14: Property landscape of the 80k imines that are in applicability domain of the consensus model. The map resolution is 31x31 and the number of RBF is 19x19. Each node is colored according to the mean logKeq value of all the compounds that reside in it. Red lines delineate the zones where "reactive" compounds are located, while the blue line shows the zone populated by "inert".

Figure 6 . 15 :

 615 Figure 6.15: Examples of "inert"-amines (left) and "inert"-aldehydes (right). Their interactions with any other aldehydes and amines, respectively, in 90% of cases lead to negatively predicted logK.

Figure 6 . 16 :

 616 Figure 6.16: Examples of "reactive" amines (left) and "reactive" aldehydes (right). Their interactions with other aldehydes and amines, respectively, in 80% of the cases lead to logK > 4.

Figure 6 . 17 :

 617 Figure 6.17: Examples of "versatile" reactants. These reactants have been found to yield imines with a very spread range of logK (from -8 to 8).

Figure 6 . 19 :

 619 Figure 6.19: Structures of aldehydes, amines and imines forming the DCL.

Figure 6 . 21 :

 621 Figure 6.21: Distribution of predicted pKi values of human CA II for 80400 imines.

Figure 6 . 22 :

 622 Figure 6.22: Speciation of a hypothetical DCL of 2 aldehydes and 2 amines. The concentrations of resulting imines before/after the addition of human CA II are shown.

Table 4 -

 4

			DUD dataset	ChEMBL dataset
	ChEMBL ID	Target Name	Activ e	Inactive Active Inactive
	1827	Phosphodiesterase 5A	170	25334	691	1515
	1952	Thymidylate synthase	63	6113	124	455
	251	Adenosine A2a receptor	79	28001	1303	3618
	260	MAP kinase p38 alpha	100	32925	1453	2567
	279	Vascular endothelial growth factor receptor 2	94	22595	2047	4663
	301	Cyclin-dependent kinase 2	189	25675	638	2305
	4282	Serine/threonine-protein kinase AKT	52	14228	725	2619
	4338	Purine nucleoside phosphorylase	102	6334	100	111
	4439	TGF-beta receptor type I	82	8013	282	385

1:Description of target-specific subsets used for model training (ChEMBL) and VS (DUD).

Table 5 - 2 :

 52 29 confirmed hit structures, experimental values of ∆Tm and IC50 (if available)

	Structure	ΔTm(pos) IC50 (μM)
		1.1	107.41
		0.76	63.28
		2.66	11.72
		1.82	160.03
		1.90	23.87
		1.60	33.16

In silico speciation assessment of Dynamic

Table 6 - 1 :

 61 Used fragmentation schemes of ISIDA descriptors and their meaning.

Table 6 - 2 :

 62 Diversity of selected libraries in 5 used descriptor spaces. The diversity has been measured using "All-Soergel" and "Min-Soergel" scores. The first calculates the average of the pairwise distances for all the compounds; the latter calculates the average of the distance from each compound to its closest neighbor.

			Product-based	Reactant-based	Cell-based	Random
			All-	Min-	All-	Min-	All-	Min-	All-	Min-
			Soergel	Soergel	Soergel	Soergel	Soergel	Soergel	Soergel	Soergel
		Mean	0,874	0,428	0,886	0,179	0,628	0,158	0,497	0,118
	IA-FF-FC-2-3	Std. Dev Median 0,874 0,002 Min 0,869	0,003 0,428 0,419	0,005 0,889 0,874	0,003 0,179 0,172	0,011 0,629 0,591	0,006 0,158 0,141	0,024 0,495 0,463	0,010 0,119 0,106
		Max	0,881	0,435	0,893	0,187	0,655	0,173	0,547	0,139
		Mean	0,767	0,378	0,812	0,132	0,628	0,142	0,569	0,135
	IAB-FC-1-6	Std. Dev Median 0,767 0,003 Min 0,761	0,002 0,378 0,373	0,004 0,810 0,807	0,005 0,130 0,125	0,010 0,604 0,604	0,005 0,131 0,131	0,028 0,583 0,506	0,006 0,137 0,122
		Max	0,774	0,383	0,821	0,139	0,651	0,153	0,604	0,141
		Mean	0,787	0,399	0,825	0,173	0,610	0,158	0,525	0,150
	IIRAB-FF-1-3	Std. Dev Median 0,788 0,003 Min 0,784	0,002 0,399 0,394	0,005 0,827 0,811	0,006 0,172 0,161	0,013 0,609 0,587	0,006 0,158 0,143	0,013 0,520 0,507	0,006 0,152 0,138
		Max	0,792	0,402	0,831	0,183	0,644	0,171	0,547	0,156
		Mean	0,730	0,312	0,790	0,128	0,568	0,122	0,480	0,096
	IAB-FC-2-4	Std. Dev Median 0,730 0,003 Min 0,721	0,002 0,312 0,307	0,005 0,792 0,776	0,006 0,128 0,119	0,011 0,569 0,538	0,005 0,123 0,109	0,021 0,481 0,455	0,007 0,097 0,087
		Max	0,736	0,316	0,794	0,141	0,591	0,138	0,521	0,107
		Mean	0,775	0,354	0,814	0,172	0,504	0,112	0,467	0,106
	IIRA-FF-FC-1-2	Std. Dev Median 0,775 0,003 Min 0,766	0,002 0,354 0,349	0,006 0,813 0,799	0,004 0,173 0,162	0,017 0,505 0,446	0,006 0,113 0,095	0,021 0,470 0,424	0,003 0,106 0,101
		Max	0,780	0,359	0,824	0,177	0,541	0,128	0,497	0,111

Table 6 -

 6 Data coverage provided by a diverse library. In the context of this study, the chemical space coverage is defined according to the fragment control applicability domain approach.

	Figure 6.7:
			Product-based Reactant-based	Cell-based	Random
		Mean	92,48	32,66	98,78	99,53
	IA-FF-FC-2-3	Std. Dev Median Min	0,65 92,21 92,21	6,13 34,47 20,81	1,46 99,17 95,85	0,81 100,00 97,92
		Max	94,24	39,02	100,00	100,00
		Mean	98,81	33,63	99,72	98,67
	IAB-FC-1-6	Std. Dev Median Min	0,34 98,92 97,83	5,68 34,13 24,27	0,49 100,00 98,63	1,46 99,28 95,93
		Max	98,92	38,92	100,00	100,00
	IIRAB-FF-1-3	Mean Std. Dev Median Min	93,06 2,24 93,18 88,18	37,68 3,78 37,31 36,59	98,05 1,71 98,99 95,36	98,31 1,05 98,92 96,35

3: Data coverage (%) of the initial dataset of compounds provided by selected diverse libraries.

Table 6 -

 6 4: 5-fold cross-validation performance of the individual logK (in chloroform) regression models that form a consensus model. Fragmentation scheme nomenclature in column 1 denotes the fragment type (I-sequence, II-circular fragments), the nature of captured information (A-atom types are captured, B -bond orders are captured), the coloring scheme (FF -force field type-based labeling supersedes default labeling by atomic symbol), other options (FC -formal charges are considered). Experimental vs. predicted logK values plot of consensus model obtained from 16 SVR models (see Table1-4). The model's performance is R²=0.93 and the RMSE=0.62 log units. A gray dotted line corresponds to ideal predictions.

	Descriptor space	R²	RMSE	Data coverage of the initial set of 120k compounds (%)
	IAB-1-3	0.92	0.66	66.56
	IIAB-1-2	0.92	0.65	52.19
	IA-FF-FC-1-2	0.92	0.64	32.62
	IAB-1-4	0.92	0.63	23.52
	IA-FF-FC-1-3	0.92	0.64	9.91
	IAB-1-5	0.92	0.65	7.42
	IIA-FF-FC-1-2	0.92	0.63	7.34
	IA-FF-FC-1-4	0.92	0.65	4.18
	IAB-1-6	0.92	0.65	4.00
	IAB-1-7	0.92	0.64	2.81
	IA-FF-FC-1-5	0.92	0.64	1.69
	IA-FF-FC-1-6	0.92	0.66	1.00
	IA-FF-FC-1-7	0.92	0.65	0.93
	IIA-FF-FC-1-3	0.92	0.64	0.43
	IIAB-1-4	0.92	0.64	0.28
	IIA-FF-FC-1-4	0.92	0.65	0.27

Table 6 - 5 :

 65 The free energies of the compounds used in DCL modeling.

	Compound	Water (hartree)	Non-polar solvent (hartree)	∆G (hartree) solvation ∆G (kJ/mol) solvation
	A1	-590,293 -590,291 -0,00162	-4,30
	A2	-495,866 -495,865 -0,00154	-4,08
	B1	-557,785 -557,783 -0,00197	-5,23
	B2	-173,196 -173,195 -0,00109	-2,88

Table 6 - 6 :

 66 Cross-validation performance of the five individual pKi prediction models that form the consensus model used in this work to estimate the affinity of imines for the active site of the human carbonic anhydrase II protein.

	Descriptor space [a]	RMSE	Q 2
	IIAB-FF-1-2	0.37	0.932
	IIA--P-FC-1-5	0.27	0.963
	IA-FF-P-FC-2-7	0.47	0.889
	IIA--1-3	0.28	0.960
	IAB-FF-P-2-6	0.40	0.918

Table 6 - 9 :

 69 Speciation of a DCL formed of 2 aldehydes. 2 amines in the presence of an effector as a function of logK of the four imines. The pKi values and the proportions of free (complexated)

	imines correspond to the imines in following order A1B1/A1B2/A2B1/A2B2.
	logK of 4 imines	Proportions of free (complexed) imines (%)
	4.25/4.25/4.25/4.25	5.4 (70.8) / 17.9 (4.1) / 17.9 (4.1) / 58.7 (13.6)
	4.0/4.25/4.25/4.25	4.6 (66.9) / 21.0 (5.5) / 21.0 (5.5) / 54.0 (14.0)
	3.75/4.25/4.25/4.25	3.8 (62.3) / 24.3 (7.1) / 24.3 (7.1) / 49.0 (14.3)
	3.25/4.25/4.25/4.25	2.5 (51.4) / 31.0 (11.4) / 31.0 (11.4) / 38.7 (14.2)
	3.0/4.25/4.25/4.25	2.0 (45.3) / 34.3 (14.0) / 34.3 (14.0) / 33.6 (13.8)
	3.95/4.25/4.25/3.95	3.6 (61.1) / 25.6 (7.8) / 25.6 (7.8) / 45.7 (13.8)
	3.65/4.25/4.25/3.65	2.2 (48.8) / 33.6 (13.2) / 33.6 (13.2) / 32.2 (12.6)
	4.0/4.5/4.25/4.25	3.9 (62.6) / 24.8 (7.2) / 24.2 (7.0) / 49.1 (14.2)
	3.75/4.75/4.25/4.25	2.6 (52.2) / 32.2 (11.6) / 30.1 (11.1) / 38.6 (13.8)
	3.5/5.0/4.25/4.25	1.7 (40.5) / 39.4 (17.1) / 37.6 (16.4) / 28.3 (12.3)
	6.2.3.

2 Speciation simulation in model DCL based on predicted logK and pKi values.

Table 6 - 10 :

 610 Predicted logK and pKi of the species present in DCL.

	Compound	logK in water	pKi
	Aldehyde A1	-	4.63
	Aldehyde A2	-	5.57
	Amine B1	-	6.13
	Amine B2	-	3.54
	Imine A1B1	4.17	6.11
	Imine A1B2	4.26	4.63
	Imine A2B1	3.84	6.06
	Imine A2B2	3.96	4.99
	The modeling of the dynamic behavior of this DCL shows (

Table 6 -11, Figure 6.22) that

 6

before the addition of the human CA II the concentrations of all the imines in solution are almost equal. According to the predicted pKi values, the imine A1B1 has the highest binding affinity with the effector (pKi = 6.11), but binding with its antagonist A2B1 is competitive (pKi = 6.06). Initially, the concentrations of all the species (reactants + effector) are set to 10 mmol. Following the Le Chatelier principle, it is expected that the effector would mostly bind A1B1 and A2B1, thus leading to a shift of equilibria in DCL toward these imines. In turn, this would decrease the concentration of free reactants A1 and B1 in solution. The calculations well reproduce these effects, see Table 6-11.

Table 6 -11:

 6 Speciation of the DCL before and after the addition of human CA II to the solution.

	Note that initial concentrations of the two aldehydes, two amines and the effector were set to 10
	mmol/L.
	Compound	Concentration (mmol) Before After
	Aldehyde A1	0.44	0.34
	Aldehyde A2	0.88	0.64
	Amine B1	0.73	0.21
	Amine B2	0.60	0.66
	Imine A1B1	4.77	1.08

References ............................................................................................. 157

Figure 4.2: Activity landscapes for ChEMBL1952.

Figure 4.3: Activity landscapes for ChEMBL251.

Figure 4.4: Activity landscapes for ChEMBL260.

Figure 4.5: Activity landscapes for ChEMBL279.

Figure 4.7: Activity landscapes for ChEMBL4282.

Figure 4.8: Activity landscapes for ChEMBL4338.

Figure 4.9: Activity landscapes for ChEMBL4439.

Acknowledgements

The information on the input compounds (molecule ID and molecule name)

and the results like the consensus prediction (Predicted value field), number of applied models and the prediction confidence will be shown on the new webpage, if needed, they can be downloaded.

= 100 ⁄ ⁄

Where Actives100 is the number of true positives in the top 100 compounds, Activestotal is the total number of active compounds in the dataset, Ntotal is the total number of compounds in the target-specific dataset.

Detailed description of the results is given in the article published in J. Chemical

Informatics and Modeling, see below.

Supporting information

Supporting information includes a list of targets used for building of the universal maps (Table 4-2) and activity landscapes for 9 selected ChEMBL targets (Figures 4.1 -4.9).

Activity landscapes for nine studied subsets of ChEMBL targets selected from the DUD database are presented below. Red zones are exclusively populated by active molecules, blue populated by inactive molecules, whereas yellow and green colors characterize zones populated by both active and inactive compounds.

Supporting information

Supporting Information includes information about (i) developed pharmacophore models (

100

A typical example of DCL is imine formation from aldehyde and amine (Figure 6.

2).

Since this reaction is reversible, mixing m aldehydes with n amines resulting in the formation of mxn imines with different combinations of R1 and R2.

In DCL formed by two aldehydes (A1 and A2) and two amines (B1 and B2), four products are expected: A1B1, A2B1, A1B2 and A2B2. As one may see from Figure 6.3, two pairs of "opposite products," A1B1/A2B2 and A2B1/A1B2, are in an agonistic relationship (green lines), whereas "adjacent products" are in an antagonistic relationship (red lines). The products A1B1/A2B2 in the agonistic relationship favor the formation of acknowledged being potentially useful for biological activity predictions. Also, while the large set of affinity data supported the default "aggressive" 12x repeated 3-fold crossvalidation scheme prone by the used GA optimizer tool, logK data is less robust and was subjected to 5-fold cross-validation.

Modeling of equilibrium constants of imines formation

Imine formation data

The training set for imines formation was selected using MaxMin algorithms applied to reactants subsets, as explained in section 6. Experiments with no detected measuring issues were labeled "Exact". "Estimated" label has been assigned because of (i) too weak concentration of reactants/products or (ii) peaks superposition, which leads to the difficulties in quantitative identification of compounds.

Modeling of logK in water

Since interactions imine -biological target (effector) proceed in aqueous solution, we need to estimate logK in water using the values predicted or measured in chloroform. For this purpose, a thermodynamic cycle shown in show that the equilibrium concentrations of the products are no longer quasi-equal if :

• logK of only one imine changes on > 0.25 compared to its initial value,

• logK of agonistic imines either both increase or both decrease on > 0.15

• logK of antagonistic imines increases for one and decreases for another one on > 0.10

logK of 4 imines Proportions of obtained imines (%)

Or in the case of an input file: