Scientist Stories: Jean-Philippe Vert, Machine Learning for Single Cell Omics

Good morning, everyone. This is a new data science seminar and Alice Live joint event today. Welcome to our speaker, Jean-Philippe Bert, from Mines ParisTech and Google Brain Paris. And also a warm welcome to Julio Sáez Rodríguez, who will introduce Jean-Philippe Bert.

Thanks, Doreen. So indeed, it is my great pleasure to introduce our speaker today, Jean-Philippe Bert. As Doreen said, he's a research scientist at Google Brain in Paris as well as adjunct research professor at Mines ParisTech Center for Computational Biology, and he's also a fellow of ELLIS. So Jean-Philippe holds a PhD in mathematics from Paris 6 University, and then after doctoral training Kyoto University, he was back to Paris where he was founding director of the Center for Computational Biology, and then a team leader at KubeInstitut. And he also holds distinct positions in Berkeley and was also professor at École Normale Supérieure in Paris. And as we will find out now, he has a strong interest on the application of machine learning methods to high-dimensional and structural data, and in particular to single-cell omics data. And six months is titled, which applications that also go into disease and in particular to look at cancer heterogeneity and response to therapy. So without further ado, Jean-Philippe, thanks so much for meeting us today and we look forward to your talk.

Thank you. Thanks so much for the kind introduction and the invitation to this seminar. It's a great pleasure. I wish I could come to Heidelberg, but maybe next time. So yeah, I'm gonna— by the way, how long, how long is the seminar? Do we plan for 1 hour including questions? Yes, 45 minutes and 15 minutes for the question and answer section. Okay, it works. Thanks. All right, so today I will present some work I did and also refer some other works related to use of machine learning techniques for the analysis of single-cell genomics data. And so maybe to introduce or motivate this topic, I mean, I guess here most of you are aware of this, but there's been huge progress done in using what we call bulk genomics, so bulk referring to the fact that we can do genomic analysis including DNA sequencing or transcriptomic analysis, etc., at the level of tissues or samples. And over the last two decades or so, we, you know, we have used that a lot, for example, in cancer research to identify subpopulations of cancer or these kind of things. So this is great. However, there are cases, in particular in developmental biology or in cancer research, where we're not only interested by, you know, one genomic profile for a tumor, but we're also interested to decipher the complexity of the tumor. Oh, I see some messages that my voice quality is bad. Let me check quickly, I'm sorry for that, if there is anything I can do with that. Okay, so I'm using a different mic. Which is better? Thanks. All right, so much better. Yes, yeah, so sorry for that. So there are cases, as I said, where instead of having just one profile for a tumor, we want, for example, to describe the complexity of the tumor, which may be made of many different types of cells, even though they are all, you know, human cells, but you may have different types of cells like immune cells, cancer cells, normal cells, fibroblasts, etc. And even within cancer cells, many of you may have heard of tumor heterogeneity, meaning that even within cancer cells, we're interested sometimes to decipher how subpopulations of cancer cells— because the tumor is an evolving organization. So how do we manage to decipher and to explore this complexity within a tissue? Well, one good candidate for that is to replace the bulk genomics by what's called single-cell genomics, which basically uses kind of the same techniques but with a different protocol to prepare the sample. Long story short, which allows to profile, so for example to do DNA sequencing or to do RNA transcriptomics, not over millions of cells aggregated together, but over individual cells within a sample. So in particular, instead of— when you have— I'm more interested myself in cancer research— when you take a cancer biopsy or a small sample of cancer, then instead of computing one profile to summarize the expression of genes in the tumor, now we have a matrix of data which tells us— well, imagine that In the matrix, each column of this matrix would be one cell in the tumor, and so each column would tell us how the genes are expressed within one cell of the tumor, and this for many cells within the tumor. So this means that we are able to maybe detect that not all the cells are equal, and to do some analysis on this data in order, for example, to do things like the ones on the right. On this slide that show an image that tries to illustrate the fact that there is a big diversity of cells, and then that could be the start of some further analysis of subpopulation of cells, for example. Okay, another interesting thing about single-cell genomics is that the throughput of the techniques and the, let's say, the number of cells which can be analyzed in parallel, as well as the quantity of information we collect for each cell, It's quickly getting to big numbers. So in this— I'm sorry, I don't cite the paper— this is taken from a publication a few years back that shows how many cells we can do in one study. And you see, I don't have more recent results, but it only goes up to 2017, and already there in just a few years, the technologies have improved so much that it was possible to assay in parallel tens of thousands or hundreds of thousands of cells. I think now it's quite common to do hundreds of thousands or millions of cells in parallel in these studies. Probably it's possible to do more, even though maybe the usefulness of doing more sometimes could be questioned. But you see that suddenly we're able— if we talk about this matrix, it means that one experiment can provide matrices with hundreds of thousands of columns, which are the cells, and many rows that could be, for example, all the genes or all positions in DNA if you do genomics, for example. So these are big numbers, and for statistical methods, this is good news, because if you want to do statistics or machine learning, let's say data-driven analysis, it's usually helpful to have more data than to have less data. So since we're entering in the domain where, you know, good techniques for large data high-dimensional large datasets could be useful to analyze these biological samples. All right, so this is a bit general, but now I will deep dive a bit into 3 topics which are, you know, arbitrarily chosen because these are topics on which I did some work with some collaborators. It doesn't mean this is all that can be done or that people do on single cells. This is just what I did, to illustrate a bit the techniques or the questions the techniques were developed to come up with answers. So first topic would be about how do you move from the raw data that the machine produces to some representation that hopefully captures biology and can be used for downstream tasks. Second topic will be a bit specific about how can you infer gene regulatory networks from single-cell data. And the last topic is an ongoing work about how do you integrate multi-omics data, meaning when you do single-cell genomics over different modalities like DNA sequencing and gene expression. And so I'll start with the first part, which is extracting signal from raw data. And so if I come back to this picture that describes the protocol, so saying you put your sample, I mean prepare your sample, you do sequencing, you obtain some raw matrix, and then you do some analysis. We focus on this part, which suggests, you know, that there is something needed to go from, let's say, a matrix of numbers that would tell you, well, in— so each column would be one cell, each row, think of it as a count, like how many DNA reads have you seen for a particular gene in a particular cell. And then the question becomes, how do you transform that into some representation? So in this case,— this visual representation, one cell is just one point in two dimensions, but it could be a bit more complicated. But the goal is that the representation extracts the information from the raw count matrix in order to then be able to do something useful like identification of populations of cells. Now the reason why it's not completely trivial is related to the nature of the data. And so for those of you I mean, I suspect that for many people in the audience this is quite well known, but for those of you who have never seen or manipulated single-cell data, sometimes it comes as a shock when you see them the first time. For example, this is a snapshot of real data, which is a gene expression measure of single cells. So in this, you know, this is the raw data. So in each— on this table, each column is one cell with some arbitrary name. Each row is a gene, and the numbers represent how many DNA reads we have sequenced coming from a particular gene. So it's supposed to be proportional to gene expression. And the first striking phenomenon that you see is that it's mostly zeros, right? There are a few non-zeros, but most of the values are zero. So this is something well known, the classical characteristic of single-cell data. So in this case, gene expression data. I think often we have like 80 to 90 or more percent of zeros in the table. And something important to have in mind is that these zeros are not always— so do not always mean that the gene is not expressed, but often they just mean that the gene was not captured by the protocol. And one reason, one important reason maybe to just understand that is that, you know, doing more cells is good. So when you do single cell, you have access to the signal on any cells, but somehow the total number of counts in one experiment is fixed. Like the total number of counts in your matrix is a constant because all the counts correspond to reads which are sequenced in a machine. And so if you go to more cells, you have to have less reads per cell. So in a sense, there is a trade-off between how many cells you can assay in parallel, and what's the amount of information you can capture in each cell. And so, you know, this amount of zero is one way, I mean, one consequence of that, which is that in order to have many cells, the price you pay is that you can have less information per cell, plus the fact, you know, that within each cell, I mean, we are talking of a few molecules, so it's hard to capture them. So the price you pay is that you have less signal in each individual. So this means that, you know, if you just did a single-cell experiment like this and then you take one column, say I will study this particular cell, then it's not a very good idea because the signal in each column is very weak, very— it's mostly zeros. And so that you need to do some joint analysis of the full matrix in order to correct or infer the missing values, for example. All right, so number of zeros is one aspect, but there are a few other aspects which were documented already a few years ago. For example, it's well known that not only you have many zeros, but when you have non-zero values, then you have what's called overdispersion, you have outliers, you have strong batch effects. Like, if you do twice the same experiments on two machines or with two different operators, it's very hard to get similar matrices. So it's not that, you know, to summarize, these data are super interesting, but but they are a bit hard to believe. I mean, you cannot just take the raw data and just say, I would do a clustering directly on my raw data, because you would just capture mostly technical effects. You need to be careful. So it turns out that the question of going from the raw data to some representation— for example, when I talk of representation, the most standout one I think so far is to say, let's try to represent each cell by lower-dimensional vector, maybe not 2-dimensional but maybe 50-dimensional? How do you go from the count matrix to a smaller matrix where you have as many columns, less dimensions, and, you know, information in the dimensions? So this operation has been the topic of much research in the recent years. Here I just mentioned— I'm sorry, I probably missed many important references, and I'm sorry if some Some persons in the audience have developed methods which are not here. The goal is just here to say that, well, first, all of these papers are quite recent. Two, they use different technologies, and what you see here just in the titles is that many of the methods use machine learning-based approaches. Here we see single cell, so we see deep count autoencoders, variational autoencoders, Poincaré maps, etc. So I don't know if it's just because the community is going in that direction, or because it's so useful for the data. We try to provide some answers to that, but you see there's a surge of new methods which are developed in order to process this data with lots of zeros, with real physical properties, and do something useful with them. So to go a bit deeper, let me present one of these methods. Actually, probably it's not the best one, I don't think it's the best one, but at least that's one I co-developed, so I'm I know it well, it's called ZinniWave, and my goal here is just to illustrate the idea behind this method, and you will see that it's common to other methods. So the idea here is that you have some counts data, think of your count matrix, and you have a lot of particular statistical properties, bunch of zeros, overdispersion, etc., and you want to extract from that the signal, like imagine that the the cells live in a low-dimensional subspace, so you want to infer the position of each cell in this subspace. And so in this model, CinduWave, which is by the way some project I did when I visited Sandrine Dudroit at UC Berkeley with fantastic collaborators David Erriso, Fanny Perrotteau, and Zvetlana Grigkova. So in this model, what we do is quite standard, I would say. It's It's a statistical model where we say, well, we have count data, so we see this count matrix as the realization of a random variable, basically, and we try to model the law of the random variable. So, for example, this ugly equation here tells you that in order to observe a count y, then there are two components: one which is the probability that it's a zero or not, because we know there are some zeros that you observe which are zero just because you missed the capture of the DNA, of the RNA, for example, and if it's non-zero, it's count data, so you need to describe a distribution over integers. Natural candidate is a Poisson model, but because there is overdispersion, we move to what's called the negative binomial, etc. I don't want to lose you in the details, but just want to insist on the fact that there is lots of modeling in this part, so how do you model the data you observe, trying to put some information about what we've observed about over-dispersion, over-dropouts, over-excess of zeros, etc. And then what we need is to relate this distribution to some parameters, and this is the goal of this image here, that tell, you know, that we model in our case, like we need to model the probability of zeros, we need to model the expectation of counts when there are not zero in a given cell for a given gene, and this is where we look for a low-dimensional representation. It's a parametric statistical model where the parameters correspond to the representation of cells we want to infer, and then we just fit this model to the data by maximum likelihood in order to infer this low-dimensional representation from the count data. But so you see, it's a bit more complicated than, say, a standard PCA directly on the count matrix. Here we go through this modeling step steps in order to say that there is a link between the low-dimensional representation of the data and the observed data through this ugly link function, which can be implemented and requires some optimization. So is it useful to do that? Well, on that paper, which again is not the most recent one in this family of papers, we had a few stories explaining why it's beneficial to do all that. In particular, I would just mention on this slide, the bottom right plot, which is an analysis where, you know, we say, well, suppose you— so we take real data where we know there are different subpopulations of cells, and we check if after embedding the cells in low dimension, we can discover subpopulation of cells just by doing clustering. And we do that when we increase when we have more and more zeros in the matrix. So this is what happens if you want to have more cells, for example. So long story short, many methods like PCA or the method at the time work well when you don't have too many zeros, which is on the left of these plots, like 40% of zeros. But when you move to 80%, which is quite common in practice, many methods do not work as well. Like, there is really a cost. It's very hard when you have more zeros. You have less signal,, and it's harder to find a representation where the subtypes of cells are present. And in this case, we compare the few basic methods to the Zinn-Biwa model, which is the only curve on the top. You see that Zinn-Biwa is more robust to having multiple classes. So I don't want to oversell the method, and in particular I think there are other methods now which are quite competitive, but just wanted to highlight that the reason, for example, in this case for the robustness to increasing zeros is probably because we explicitly model this in the model, and we infer the model taking into account that we want to be able to generate zeros, so the model is not perturbed by the presence of zeros. It models it, and it's, for example, very different from PCA that would see zeros and try to, you know, find a representation to explain the zeros directly using an equation. So, as I said, you know, there are many other papers published after the EasyB-Wave. Maybe I can cite one of them, and just to highlight that it's more complicated. So this is a great work called STI by Roman Lopez at Berkeley, and Mike Jordan and Yaa Youssef. So I will not go to details or equations, plots, etc., but the only thing I will mention is that even though the architecture is more complex, it uses deep neural networks to map data to some low-dimensional thing, that can then be back to data. What's interesting again is that if you look at, you know, besides the architecture, what's the objective function, it's a bit similar because what it does is that from the low-dimensional representation, it learns a function that learns different parameters which end up being expected counts, expected dropout, etc., and then models the count by You know, it's written in a different way, but it's basically similar to the zero-inflated negative binomial that I presented in ZimB. So this means, to summarize a bit, that there is a commonality to say we want to model data with some specific plot, model the particularities of the data, and then different teams have proposed different models. Here we just work on the low-dimensional factorization, matrix factorization. SCVI uses a different mechanism, which is deep learning to learn a representation, but at the end of the day, in both cases, the goal is to model some distribution of the data that then can be fit by maximum likelihood. Okay, so because there is so much, one question that users can have is: should you use PCA, should you use CBWave, should you use SCVI, or something else? That's the question we try to address with a student of mine, Félix Raimundo, in collaboration with Céline Vallot from the CURE Institute. And so we did the benchmark, which we published last year, on just the question of— we did not benchmark all methods, but 5 methods, including source-gram is based on PCA, is easily represented, SCVI is the one I mentioned, there are a few others.. And one question was, well, do they work differently? Do they provide different representations? And if yes, which one is the best one? And I said, so we did the benchmark on 10 datasets where we were able to play a bit with the mix of cells, and we focused principally on when you do this embedding, can you then recover populations of cells? If you know, for example, that you mix two different populations of cells, can you, through these single-cell experiments, recover them easily. And so I will not describe all the results of the benchmark because, you know, not too many specific results, but maybe the main one is that it seems that if you take the methods with their default parameters, which is what most people would do, like you download the package, you run it, you obtain results, then In many cases, in particular in the hard cases, which are summarized— it's hard to read, but on the plot on the right here, each set of bars is one dataset, and vertically some performance, which is how well you would cover the subtypes. You see that there are a few hard cases, so we designed a few hard benchmarks which are a mixture of 5 populations of cells, for example, or 8 populations of cells.. And if you take the method with default parameters, it's a bit surprising how good simple methods like PCA are compared to more complex methods, even like SymbiWave or SCVI, which are not very good with default values. However, if you allow to tune the parameters, because, you know, the deep learning or SymbiWave have a few parameters that you can tune, suddenly you can get better results. So I think it's quite coherent that, you know, given that more complex methods can be more powerful if you tune them well, then you can get slightly better results. But in particular, so in these datasets, which are not huge— I think the situation would change once we move to datasets with millions of cells, but here these are standard datasets with a few thousand cells— then it seems hard to beat PCA. So one conclusion of the paper is that until we have bigger datasets, maybe it's safe not to— I mean, maybe it's safe to use all the methods existing, but to be careful that sometimes the fancy methods, especially if you run them without trying to optimize the parameters in some way, can lead to disappointing results. So to me, I mean, this is just a signal that No, there is some research going on, and the community as a whole probably has to be a bit careful about, you know, what methods are used, and for example, not being too excited by fancy new methods which promise a lot, but which at the end of the day you need to be a bit careful when you use them. Okay, I'll stop it on that topic and move on to the second topic. Which is completely different. It's a bit specific as well, but it's the question of when you have single-cell data, can you use them not to learn a representation, but then to do something a bit more ambitious, which is rediscover how cells work. In particular, rediscover regulatory networks, gene regulatory networks. And so this is a topic, you know, maybe Rubio could give a much better talk than me on. That topic as a whole. How you infer regulation from expression is a broad topic, and people have worked on that on bulk data a lot. Something interesting is that when you have single-cell data, you can imagine that because you— suppose you have a dynamical process, that could be cell cycle, or that could be differentiation, or that could be tumorigenesis, these kind of things— then potentially because your sample I mean, if you are able to sample cells at different positions, let's say in time and in the differentiation process, then you may be able, just like on this picture, to infer that there is some trajectory going on. And if you have trajectories, if you are able, for example, to order the cells, then if you look at their gene expression, you may be able to fit a kind of dynamical model on the expression that then you can associate to a gene regulatory network. For example, if you observe that at some position a gene is going up and this is always the case when some other gene is active, then you may infer that the active genes may be a regulator of the gene that goes up. So the single cell provides you with these— I mean, potentially single cells can be very useful for that because they give you individual cells, and therefore if you are able to reconstruct the trajectories, then you may be able to we can constrain the dynamical equations as well. So of course, to do that, you know, if you want to do that in practice, you need first to have data where there is some dynamical process going on. Let's talk about differentiation, for example. And then you need to be able, when you do your single-cell experiments, you know, to associate each cell with the time. So I will not discuss this question, but suffice it to say that this is something that people have worked a lot on. Like if you have a population and you know there is a chemical process, how do you reorder cells and try to, you know, control what's called the pseudotime, so to give them some biological time in the process. I think it works more or less, but we assume that it works well. A bit I summarize in this paper by Paul Trappnell and colleagues, that, you know, where you have for example two experiments, each of them with a mix of cells, and it's possible to not to take the real time of the experiment, but to reinfer some ranking, sorting of the cells from the less differentiated to the most differentiated, and then you can— could be the basis then of some dynamical model. So if we assume that this is done, I would focus on the next question, which is once you have ordered your cells and you assume that there is a dynamical process going on, how can you infer the dynamics, and in particular Could you then reconstruct some gene regulatory network if you made a very strong— and there are many questions about that— the very strong assumption that just looking at the profiles of gene expression, you can infer the dynamics, so that tells you how the expression varies as a function of how the expression is, and that this could be related to regulatory network. Of course, it's a strong assumption because we forget about the link between gene expression and protein level, then all the protein regulation, etc., is not seen through the lenses of transcription. But we made this assumption to see how far we can go in just from the expression because gene regulation— So if, again, assuming that we have been able to assign a pseudotime to each cell, that the most natural and simplest dynamical model, maybe a bit naive, is a simple ordinary differential equation that tells you that at time t, if you are able to observe the vector of expression, so x would be, you know, your transcriptome profile of a cell at time t, then we would assume that from this profile it's possible to measure the dynamics as well, if we assume again that the expression of a transcription factor indicates how active the transcription factor is, and therefore how it will influence the growth or the reduction of expression of its target. So the simplest model would be to say that the derivative in time, which I write as x'(t), so let's call that the velocity, the velocity of a gene, like is it going up or down, depends on the expression of its transcription factors, matrix A. And so A, you have to think of it as a big matrix of size number of genes by number of genes, which is very sparse because basically it just says you that to measure the velocity of the, you know, how one given gene changes only depends on the expression of a few of its regulators. So the structure of the matrix A will contain the architecture of the regular tree network, and the goal now is from observations of expression, how to infer matrix A. And so assuming it's sparse, so if it's sparse, then we can from the matrix A infer a regulatory network that would tell you the non-zero entries in A, would tell you which transcription factors regulate which target genes. Okay, so one difficulty, and if I come back to this thing, you know, one difficulty is that we don't observe directly x'(t), we just observe x for many cells, so the question is how do you infer the matrix A from the observations of many points, which is not completely trivial. So one approach to that was proposed by Matsumoto and colleagues in 2017 with the method SCODE, which works as follows. So it really follows, you know, the same equation as the one I showed in previous slide, and it says, well, I don't— you know, in order to fit the matrix A, let's try to fit the matrix A that can explain to explain the data, and to explain the data, I will optimize a criterion, so find a matrix A that minimizes some cost, and the cost here is for each cell, you compare the expression that you observe, which is x(t_i), so x(t_i) would be the vector of expression of the cell at time t_i, and you compare it to what it should be if the dynamics was following the the equation here, so x prime equals Ax, starting at some position x naught. And the reason why it's possible to do that is that if, you know, the equation is so simple that if the equation is correct for a given A, then you know that if you start at the position x naught, which has to be optimized, then at time t you can mathematically solve the fact that at time t the the vector expression x_t should be exponential tA times x_0. So this is the solution, the mathematical solution of the differential equation, a very simple one. And so this compares the observed value to what it should be if the equation was correct. The question is, can you find a matrix A such that this is minimized, so that the matrix A that could explain most of the expression. Now, you know, this was a software doing that, There are a few difficulties here, which is first that this is hard to solve. So if you want to solve this problem, we don't know how to find the best matrix A because it's a non-convex problem, so there is a heuristic provided in the code. And another more conceptual problem is that here you see that, I mean, in this criterion in particular, because there is this exponential term, the solution of the ODE is very sensitive to A, in particular for large times., and so this basically, in this objective function, only the cells with large time will influence a lot the choice of a, because they can be very wrong. Like, if you make a small error on a and you solved it for a long time, then the solution of the equation at time t, when t is large, is very sensitive to the parameters of the model, which is not the case for the cells at small t. So, you know, this equation creates some some discrepancy between small time and large time, and it's not completely clear why inference of neural network, of the generative neural network, should depend mostly on the cells with large pseudotime. In particular, choice of 0 is a bit arbitrary. Okay, so with Pierre-Cyril Aubin, a student of mine, we worked on a different model which starts from the same equation, but which says, well, instead of solving the equation and then comparing the observed x with the predicted one, can we just more simply try to find a matrix A that explains the velocity at any time compared to what it should be? So, you know, x prime equals Ax was the equation, and here we just say, can we find an A that minimizes the predicted velocity to the true one? So if we are able to solve that, that would be I think, a better model because here we just focus on whether the equation— so directly fit the parameter of the differential equation to how it should relate, let's say, position in the space of transistors and velocity. It's not sensitive to large t, but the big caveat here is that in order to solve this equation, you need not only to observe x, which is the expression of a cell which should you would observe, but you need also to know the time derivative of it. So for each gene, you need to know how much it is expressed and how it is moving in time, like will it go up or down. And this is not obvious, it's not directly observed in the data you have. So in order to apply this method, you need a way to infer the velocity. But so once you say that, you say, well, maybe I can infer velocity. So for each cell, can I infer not only what's the what the expression of the genes are, but as well what are they going to help with that. So in the paper with Pierre Surer, we developed a particular method which is based on the fact that because you have many cells, you know, just from looking around, you may get a sense of what is the trajectory of the cells. And, you know, forget about the ugly equations here, but intuitively, just by looking at the cells around you, it's possible to infer basically where the cell goes, assuming that the cells form a trajectory. So this can be formulated mathematically with some estimator, but the intuition here is just to infer the velocity of a cell. I look at cells nearby, so with nearby pseudotime and probably nearby expression, and I fit some model that infers the direction of the flow of the cells, and this could give me the direction where my my cell is. Okay, so important thing is that for my cell, the vector of expression is given by experiments, and the direction which is used— so this is the x' here— the direction is inferred from the neighboring cells. It's not directly in the data. Now, after we worked on the adaptivity, there were a couple of super interesting papers that were published by different teams that proposed another approach to also infer velocity. And here it's, you know, I refer to these papers if you want to know more, but basically these are papers where in order— so they also say that in single cell sometimes you're interested not only on what is the expression in a cell, but as well as what is the velocity of the cell. So what is the time derivative if there is a process, how the cell would evolve. And so they propose methods based on super smart ideas of looking in the details, you know, at the reads, whether the reads are in introns and exons, making the ratio, and making some models to say that if you have more reads in introns than exons, etc., then you can infer that transcription is active, etc. So I will not really describe these papers, they are super nice to read, but long story short, they also infer, you know, these small arrows. So for each cell you can from the data infer not only where the cell is but as well where it goes. And each case of the arrow is quite different from the one I presented before, because before we infer the arrow just based on surrounding cells, whereas here it's more based on the raw data, so the detail of the DNA reads which are sequenced within each cell. So these are a bit like two different ways to infer the same thing, which is inference of this velocity field. Taking any of the methods potentially gives you the access— I come back a few slides earlier— so gives you access to this x'(t) for each cell, and then allows you to solve this problem. So say, let's try to find a sparse matrix A so that I can relate the expression of a cell to its velocity, which is inferred by some way. And if I can do that, then I can extract from my matrix A some candidate regular training network. Again, good thing is that this is super easy to solve and there are tons of methods to find a matrix A. This is standard least-squared error regression. Okay, so we did that, and just let me show a few results on relatively small datasets, but, you know, we need So, this is what we did. So these are two datasets, one marine data and one on human cells, and here we just compare mostly two methods. So one is on the left, SCOD, and on the right, so Grizzly, is the method that I mentioned. The performance is in area under the ROC curve, so 0.5 is you make random precision, and 1 would be that you exactly recover the regulatory network for these particular cells, which is in databases. So overall, performance is not very high. This is a hard problem, but there is a significant improvement compared between SCOD and Greasy. So the idea to change the way regression is phrased seems to be helpful in both datasets. And here, I add Tigress. So Tigress is one method, one out of many methods, that would try to do the same. But forgetting that we have single cells which form a dynamical process. So just seeing each single cell as a possible bulk gene expression data, and so TAGRES is just a method to infer a gene regulatory network from gene expression data without the hypothesis that these are single cells. And so you see that it doesn't work well, meaning that, you know, I think it's fair to say that if we have single cells, then we should use methods specific for single cells, in particular with the idea that the cells— so the individual gene expression that we observe in single cells are points in a trajectory and not some average of some— that's not any sense. Another question we asked is: do we see differences if we infer velocity either by the method I described, so where you infer the the velocity based on neighboring cells, which is Grizzly, or if you use for example SCVDO. So SCVDO was one of the methods that infers velocity from the ratio of introns and exons counts. And here we don't see big difference, but we see significant difference in favor of inferring velocity from neighboring cells. So I don't know exactly what, you know, it's not huge, but It means at least that the way you infer velocity leads to different network inference, and in this case, it's slightly better to measure velocity, to infer the velocity by the methods I described based on neighboring cells, as opposed to using the raw data. Okay, so I'll stop here with this topic, but just wanted to highlight again that It's a super fascinating topic to say, can we reconstruct the gene regulation from the raw data? It's complicated. It's not easier than it was with bulk sequencing. But we are changing the methods, still seems to be changing the performance. I think there is room for new development there. And third message is that it's quite clear that single-cell gene expression and bulk gene expression are not exactly the same when it comes to gene regulatory networks, but the methods should be exactly the same. Finally, I will finish quickly on the last topic, which is one of active research in my group and in many other groups, which is integration of multi-omics data. So here the motivation is that I talk of single- single-cell omics, and I was not very precise always in terms of how am I talking of DNA sequencing, am I talking of transcriptomics, am I talking of histone modifications, etc. In fact, nowadays with single-cell technologies, you can assess all of these modalities. However, often you cannot do that. I mean, sometimes you can, but often you cannot do it simultaneously on individual cells because, for example, once you measure on a given cell. In order to measure histone modification, you need to extract DNA, and you cannot simultaneously measure methylation, for example. You can do that on bulk, because on the bulk you split your sample in sub-bulks and you do different modalities on different bulks. But here on individual cells, once you measure one modality, often you can't measure the other one. So there's an interesting question here, which is that now that you can measure many things, and maybe it includes images as well, Can we integrate this data in order to, for example, compare what you observe in the microscope with what you observe in the transcriptome and what you observe in the epigenome, etc.? It's a hard question because, as I said, often we don't have— like here, if I take two modalities, we don't have direct— from experiments, we don't have a direct mapping between which cells are measured on gene expression and which cells are measured in ATAC-seq. Often they are not the same cells, they are different different cells, but populations are similar. So in this picture, you would like to maybe align and computationally guess that some cells in the space of RNA-seq biology correspond to some cells in the space of ATAC-seq. But again, it's hard because we don't have this from data, and also there is no— often there is no direct correspondence between what you measure in one modality and what you measure in the other one. A couple of years back, working with Bill Noble and a great team of his postdocs, we developed one method called MMDMA. I will present quickly, but long story short, this is a method that takes two matrices, so two modalities, one let's say RNA-seq, one ATAC-seq, but this could be other things, with different number of cells, different number of features, and tries to align them. But of course, aligning has to mean something, so here our definition of alignment is that we want to find a mapping of each modality in a space which would be common, so a d-dimensional space, so create these two matrices in such a way that hopefully we capture biology, and capturing biology would mean in particular that in the learned space, the two modalities overlap as populations. So the hypothesis here is that even though we measure different things on different cells, we assume that the populations of cells are the same, and therefore that there is a shared space that can be inferred from each modality where the two populations overlap. And if we are able to do that, you see that if this works, then after, at the end of the day, we have an embedding that then allows us to in silico predict, for example, one modality for the other, or match a given cell from the first sample to a cell on the second sample, the most similar cell in a different modality. So I'll skip the details of how this works because I don't have much time. It's not complicated, but it's a bit technical, and basically it's implemented as an optimization problem where probably the most important The fact is that we need to use a measure of overlap between populations, and there is lots of work in that field. We use one measure called maximum discrepancy, that could be the one, but this is a measure between two sets of points, so that is small when the sets of points overlap with each other, and that is large when they do not overlap. And one important thing is that this measure is continuous, and you can make the derivative the gradient of this method, of this measure, with respect to the position, so you can optimize the position, which allows to learn end-to-end a system that learns how to map the raw data to a space such that the two populations are parallel. So, it works quite well on toy data, and, you know, the kind of things that we did, if I jump a bit, this is typically what it does. So here, I mean, this is on a very small dataset where you have two modalities. Here we look at gene expression and DNA methylation. And so this is a dataset where we use— I mean, the data were produced by a technique where it's possible to jointly measure gene expression and methylation on the cells. This is a case where we know— so we have two matrices of data, but we know in addition which cell in one modality is which on the other modality. And we did some experiments using this data as validation to check if the algorithm was able to reconstruct this matching, which we know from experiments, but to check if the algorithms could find it if we didn't give it as input. So the kind of thing it does is, you know, on the left here, imagine that you have a random mapping of all the cells into modalities, and the dotted lines correspond to matching who we know is who. So we match hard to see maybe, but a blue point to a red point. The blue points are DNA methylation, the red points are gene expression, and we put a dotted line here if it corresponds to a given cell into modality. And so basically what MM-DNA does is, here if you map to 5 dimensions, is to learn how to map it and go from the left to the right picture where you see that the mapped cells— so each cell then has 2 points, which is how how its gene expression is mapped, how its DNA methylation is mapped. And what we hope is that if it works, then we should have, you know, the lines, dotted lines, which get small because we hope that a given cell would be mapped to the same position either from the viewpoint of gene expression or from the viewpoint of methylation. So this is just a plot, but this can be quantified. Like, if you want to assess if MNDMA does a good job or not, we can do that systematically and, for example, compute the distance between the real pairs and check, for example, if the good— so the correct gene expression given methylation is closest to— so if the pair is close to each other compared to other pairs. And so we can quantify that, and I will not detail the thing, but we have a measure which is the fraction of samples closer than the true true match, which basically measures how good we are at mapping a given cell to the same position. And we're able to show that the nMDA does a good job, in particular that in this case we see that for the majority of the cells, the nearest neighbor— so if we say we have gene expression, who would be the best candidate to be the same cell in gene methylation, DNA methylation— then in more than 50% of the cases, we're able to recover it exactly. So it means that we are able to learn something from the data. Now, there's been lots of work going on in that space in the last year or so. For example, here I mentioned this— this is not my work, this is a very great work by Kevin Wu and colleagues on one method, BABEL, that here says assume that we have data from two modalities. So this is the case where we are able to we measure two modalities, can we learn to translate, for example, from one of them? And they presented a nice deep learning-based model that's able to translate. Like, if you have the RNA expression of a cell, can you predict the ATAC-seq accessibility profile in viruses? And they gave quite interesting results. So I mentioned that because more recently, again with Bill Noble and Ran Zhang and Etisap Papadzentos. We extended a bit that work, I mean, did something a bit different, but assuming here, we assume that we have access to a small collection of cells for which we are able to do co-assays, so meaning that we're able to measure two modalities, and that in addition to that, we have access to, and this is easy to get, we have access to many cells, in other experiments where we just measure one modality. The question is, can you use the co-assay to supervise a bit something and use large datasets of single assay in order to train models that could learn to align the two kinds of points? Or we focus in particular on translation. Translation meaning that if I give you, for example, an RNA-seq profile of a cell,: Can you predict its ATAC-seq profile? And so we want to generalize a bit the BABEL model, which only used co-asseay, by asking: Can it be helpful if in addition we have plenty of cells with only RNA-seq or with only ATAC-seq? So we'll not detail the model. If you have questions, I'm happy to discuss it, but these are just models based on a deep version of autoencoders where we learn to translate between the inner layer. Long story short, it seems to be working interestingly. Maybe, you know, on the left plot here, what we show is a comparison between Babel and our model when we just use the cross-modal data. So this works well, it works similarly, I said, the two methods. And as soon as we add lots of other data which are not aligned, so in order to provide some kind of semi-supervised learning framework, so a lot of data for each involved modality,, you get a significant improvement. So here, on the curve, here, sorry, FOS CPTM is exactly this measure of fraction of cells where the nearest neighbor is the correct one. So when you are zero, it's better, and we see a significant improvement compared to BABEL or PolarBear, another method using only coarse-grained data, as soon as we add extra data. I think it's quite promising in terms of combining the expensive and technically challenging data which you obtain by QSA with the plenty of data that we are getting now. And here we're talking of, for example, 800,000 cells or 100,000 cells. So that's a lot of data with a single modality, and you see that it helps a lot here. So I think it was quite promising direction. All right, I will not detail all the other results, but maybe just mention that if you're interested in that topic, there's a competition going on right now, so don't hesitate to have a look. It's a NeurIPS competition on multimodal single-cell data integration. Could be a thing you may want to have a look at if you're interested to know more. There are some data, there is a challenge like run a model to co-embed or translate, etc. And I will stop here. So, conclusion: I tried to illustrate through different approaches that I'm personally quite excited by single-cell genomics, in particular because we start having big data. Having a matrix basically instead of a vector means that the amount of data gets large. In particular, if you want to learn things about cells, It's quite common now to have hundreds of thousands or millions of cells, and so this is exactly the kind of quantities where fancy or fancier machine learning models can have some effect. There's tons of exciting perspectives. There's many challenges. I did not insist on them, but I think this deserves lots of research on trying to better understand the structure of the data, how to trade off quality and quantity. You know, I touched upon that in the gene expression, but it's a general trend. Like, if we have better algorithms that can handle noisy data, then we can allow more cells, and this in turn can allow us to investigate, for example, populations of cells in small proportions which are often hard to capture if the signal is not good. So these are challenges. There are many things that I did not talk about. I think in the future we have also to remember that cells talk to each other What I present is a bit caricatural, saying, well, we have a sample, we extract each individual cell, we forget where they come from, we just see them as a cloud of points. Probably we need also at some point in the future to try to be better understanding not only how the cells behave as populations of individual points, but what's the communication channels among them. And technically, I think one important question is, and I think not enough people are working on that right now, is that the community of people working on mixed data coming from bioinformatics and the community working on imaging, which is also going very, very fast, and the techniques of single-cell imaging have to converge at some point. And we have, you know, we tend to have now two different approaches, either based on DNA sequencing or based on imaging. That provide signals which are the same resolutions about the same objects. So combining these two things together is likely to be of great interest in the future. I'll stop here. I'm sorry I was a bit long. I'm available for questions now, but I'm happy to stay a bit longer if you want to ask questions.

Thank you. Thanks so much, Jean-Philippe, for a really great presentation. And indeed, we have time for some questions and answers. I would say to start with those that brought it in the chat. So there were, from the first part, a couple of questions related to Zimbabwe from Sarah and Dimitri. So they were asking if it's suitable only for UMI data since it uses the negative binomial distribution.

And if the answer is yes, do you need the self-inflation to model UMI data? Yeah, it's a good question. I don't have a definite answer to that. I would say that when we developed ZimmiWave, most of the data, we did not work on UMI data, but then after that, there's been several works, you know, claiming, for example, that as soon as you move from raw counts to UMI, then the need for zero inflation is less important. So I don't have a strong opinion on that. I would say that things, so ZimmiWave or SCVI, etc., tend to use a specific model that allows to model zero inflation, and I would say that if there is no zero inflation in your data, it should not hurt the model. It would just, you know, typically the model could learn that the probability of zero inflation is zero. So it's not, you know, you can use ZB-Wave or SVI on all data, and I didn't see any publication saying that if you use them on UMA, for example, it doesn't work anymore.

I think it works.

Thanks.

There is another one on the first part. I think it's related to your benchmark paper comparison, whether complex methods discuss their own performance on parameter selection.

Some may have claimed that their performance is robust to hyperparameter selection. Yeah, so if you take— I mean, you can take publications one by one. Some do, you know, I think in some publications there is some analysis of the impact. What we did is, you know, we did a benchmark for Producible. We observed that— so to be clear, there are plenty of parameters. Some have very little effect, some have a huge effect. We observed in our benchmark that when Probably the main thing we observe, as I said, is that if you take the default parameters, which are supposed to be the good ones that people use, be careful, it doesn't work in many cases. I would say that depending on the method, sometimes the few parameters that you can play with are quite obvious. It's like, do you log-transform the counts or not? And you have to log-transform, it's quite obvious. But for example, on the deep learning-based methods, These are more things like choice of the learning rate or structure of the neural network, etc. And here we observe, so I would say, for example, that in the field of machine learning, people doing computer vision or speech recognition, they know very well that— they don't understand very well why, but they know that the learning rates of the stochastic gradient descent has a huge impact on the performance of the model, and I'm just a little bit worried that people in the field of bioinformatics, maybe who have less experience with this thing, you know, just want to take the software as it is and are not aware that they should tune the parameters carefully. So I'm not sure this is discussed enough in the publications.

I agree that for sure, at least being aware that the default parameters may not be good, that's a very important take-home message. Maybe we are at the top of the hour, maybe last question and then we can stay on the line a bit longer for sure to discuss, but if people need to leave top of the hour, there was a question around moving to the gene regulatory networks, so whether there is an assumption on the network invariance along the trajectory when the cells get different nodes.

Yeah, so in the model I presented, indeed, the My assumption is that there was this equation, x' Ax, so behind that strong hypothesis is that this is, you know, the matrix A doesn't change, meaning that in biological terms, meaning that there is a regulation, a regulatory network that is constant throughout the trajectory, and this is a super strong hypothesis.

Maybe a final question from my end. And so along the generatory network, but more in general, you know, like for example, we are very interested in how we can use different sources of biological knowledge to help the methods, right? And of course, machine learning, you can do in a pure data-driven way, or you can try to use biological knowledge. And so what is your take on, on general on this? So do you think with these huge data numbers, maybe the data alone tell us what we need, or we can use knowledge, prior knowledge, in some ways to deal with the issues such as noise and so on?

Yes, I would say probably depending on what you want to do. So for example, you know, I presented the generative network, and I mentioned that, you know, the first thing you see in the results is that the performance is poor. I think For example, inferring the full regulatory network just from data in a purely data-driven way, maybe it's not a good question to ask, and at least it doesn't work well now, and I don't see any promise that it will work in the next few years because for many reasons, maybe we don't have enough data, but we don't really have the good models, we don't know how to infer, etc. So I would say that's the case where probably The good way to phrase the problem is something more like, well, we know part of the gene network, maybe we want to infer the missing parts, or we want to infer which parts are active in some contexts, etc. And then there could be dedicated methods to test hypotheses or to infer things like that. So I think having prior knowledge may be a way to phrase correctly the problems that then lead to suitable inference methods from data.

Okay, thanks. Yeah, now I mean, if we are in kind of the extra time, but we can scan the line and people want to ask to Jean-Philippe, also feel free to speak up orally.

So anyone wants to say anything else, any other comments, questions? Jean-Philippe, I have maybe a comment. Thanks a lot for the wonderful talk. So I have a question to this final multi-ohm integration part. So you use the similarity between the sort of in the final embedding as the performance measure, right? So you You know the ground truth, you map in the expression and methylation modality, and you estimate simply by distance. Do you think this might be a too strong assumption from some systems? I don't know which data you've been using for this question, but could it be that the cells don't have to be close in both modalities? The modalities in the same cell. So there might be cells which are closer in methylation space, but far in expression and vice versa. So how would you comment on that? Yeah. So yeah, to be frank, I think that, again, different ways to phrase the problem. So we decided to phrase the problem as, Ideally, a method— maybe that's not the best slide, but something like that— what we would like is a magical model, and maybe we don't have the good one, but someone will invent a good one that is able to extract something from each modality. And the something should be, you know, what we call the "learned space" here, should be, let's say, some biological information which You can think of it as universal across modalities, meaning in particular that one property of this learned space should be that a given cell has one position in that space, and from that position, it should be able to— the learned space should have the information to be able to, for example, generate the RNA-seq or generate the ATAC-seq. So you see that in this model, by definition, one cell has one position in the learned space, and so this is coherent then with our way to assess if it works or not by saying that if we're able to go from the data to the learned space, then because each cell should have a unique position, we measure how well it works by computing how far our inferred positions are from each other when we infer it from one modality or the other. Now, this being said, maybe this is not a good way to formalize the problem. Maybe, so for example, you mentioned maybe we don't want a given cell to, you know, maybe one cell should have two positions in that space. So if we do that, it means that we need to change a bit the formulation of the problem and say that what we want to learn, for example, is a mapping from modalities to a space, that has some usefulness in some other ways. I don't have a good example in mind, but I agree with you that maybe it's possible to rephrase the problem as with the property that a given cell may have different mappings in the space, and what you want in the learned space could be something else. For example, maybe to be more concrete, you know, at the end I mentioned this polar bear work, or baboon work. I insisted on the fact that we tested it on on that, but I think the main test we wanted to do is a bit like double. The goal was to translate, so to be able to say, if I do one modality on one cell, can I predict another one? And there, maybe you don't need to have a unique latent space. It turns out, for example, that we do have two latent spaces in our model, so one latent space for RNA-seq and one latent space for ATAC-seq. They are different, but they have the same dimension, so potentially we could we can map them, and what we do is we learn a translator between these two spaces. So when I say let's compare the position of the cell in one modality to the same cell in the other modality, in fact I did not give all the details. I said let's compare the position of the cell after translation to the position of a cell in the other modality. But you see that before translation, a given cell has two different positions, and that's good because that's good enough for us to solve the problem of translating, right? So that's the case where the way to measure if a mapping is good is not through the fact that the cell would be at the same position in two modalities, but rather that we have a simple way to go from one position in one modality to what it should be in the other modality.

Yeah, thanks a lot. I think so. If there is no other question, uh, I see none in the chat. So thanks again, then, Philip.

Was great to—

you could join us virtually, and as you said, hopefully next time in person.

And yeah, thanks everyone for joining, and, uh, have a good, good day. Yeah, thank you very much. Yeah, thank you for being here. Thank you, Sean Felley, for the talk. Again, thank you very much. Have a good day, everyone.