Scientist Stories: Hani Goodarzi, Deciphering the RNA Structural Code

If you want to know more about different aspects of my science that are not covered today, feel free to reach out. But the topic that we picked for today is really focusing more on how regulatory information is encoded in RNA and how can that information then be decoded computationally or experimentally. So really, I want to take a step back and start from really the very beginning. When I started to think about this problem. You know, we all know that all the cells in the body really share the same genetic blueprint, so the diversity that we see in form and function in the cellular landscape really arises from which of these genes are expressed, to what extent. Basically, spatial and temporal mapping of gene expression patterns really is a driver behind, behind the cellular diversity. Ultimately, the idea is that there is a regulatory program in place in every cell that actually controls this process. So, you know, each regulatory program can be cell type-specific, or it can be universal or broader. So, ultimately, this is really the master regulatory program that we are trying to understand and study. It can for the most part, you know, serve as a black box, but that's exactly what we want to peer into. Now, we do this, you know, you can do this in the context of, of course, healthy cells as well, like developmental processes, for example, or differentiation models, but over the years we have looked at this problem more through the lens of cancer progression. You know, the same idea applies in that space as well, that, you know, you start from a healthy cell that has an intact regulatory program that decides the gene expression patterns and gene expression landscape of the cell, and as the disease progresses, different components of this regulatory network is actually hijacked by cancer cells to achieve dysregulated gene expression patterns that are actually driving the disease. And as you go further and further away from what is normal, you basically get more and more dysregulations. And, you know, there are basically many non-genetic factors that are driving cancer, for example, especially in the later stages. And my lab specifically is more interested in this, this final step in progression of the disease where you basically get— go from a primary tumor cell to actually metastatic disease. And we're trying to understand how cancer cells spread throughout the body and recolonize tumors in other organs. And we know that through this process, again, there is a fair bit of gene expression reprogramming as well, actually one of the hallmarks of this process. And here on the right, I'm kind of visualizing that in a simple way by comparing primary and metastatic breast cancer tumors and just really listing genes that are differentially up- or downregulated to emphasize that there are thousands of genes that are differentially expressed in a traditional sense. But really, if you take a step back and look at the entirety of the gene expression landscape, it's really— the cell is basically transformed. So we are less interested in individual events in that sense, because answer to the question— the answer to the question of you know, what are the genes that are changing in the context of a complex disease like cancer is not that helpful because the answer usually is almost everything is changing depending on the resolution that you're looking at the problem. So we are really instead interested in understanding what are the underlying regulatory programs that are responsible for this type of reprogramming in the gene expression landscape that we see. Now, of course, historically, when we say gene expression, for the most part, folks have been focused on basically the steady-state mRNA levels in the cell that is, as, you know, as you probably know, is driven by both transcription and decay. So it's basically— this is steady-state levels of mRNAs that have become a proxy for gene expression in the field. Initially, measured using microarrays and more recently with RNA sequencing. But of course, gene expression control is a cascade, it's a multi-step cascade. Each of these steps are very well regulated, you know, from transcription to RNA processing to RNA localization, either to the cytoplasm or different parts of the cytoplasm or even the nucleus. And then If it's a coding RNA, you can have translation and ultimately the decay of RNA. So then this entire process is, is very well regulated, and this regulation happens, for example, really through the action of transactors that talk to DNA or RNA through cis-regulatory elements. So, you know, at the transcriptional level, for example, you have entities like transcription factors that I'm sure everyone is familiar with that bind linear elements in the promoters or enhancers, and they are the ones are driving gene expression of a given gene. And at the post-transcriptional level, for example, you have RNA-binding proteins or microRNAs that similarly bind cis elements on the mRNA molecules themselves and impact their lifespan, for example, impacting their degradation. Now, the question that we would like to answer fundamentally, as I mentioned, is if we are given basically a data set of gene expression changes, for example, in the context of disease or really any other context, from the knowledge of how genes are changing in different conditions or in different scenarios, can we figure out what are the underlying— these underlying regulatory programs that are responsible for the emergence of these patterns? So this is not a new question. It actually goes ways back. It was an important question when I started grad school, is still an important question today. But I think fundamentally there were key advances in the field that I'm trying to cover today in a few slides and then link it to the kind of work that we have been doing on RNA structure, because fundamentally the two questions are interlinked a bit. So, uh, one of the, uh, kind of, I think, greatest insights in the field emerged from this understanding that, well, if we are talking about regulatory programs, ultimately genes that share the same regulatory programs then should behave similarly. So the opposite of that becomes if you look at genes that are coexpressed, so meaning that they show similar patterns of gene expression, in their vicinity, either in the coding sequence or promoters or enhancers, they should share similar binding patterns of transactors, for example. And this was really the basis of a lot of early— and even today the basis of basically cis-regulatory element discovery in the field of genomics. And from very early days, there were tools that were developed for basically taking a number of sequences and looking for representations of motifs in them, as I'll talk about. And again, some of these tools, you know, are going strong to this day. You know, the MEME Suite is still kind of like a pillar of the community for motif discovery. So this has been a long-standing problem, and there are a lot of existing tools that for many years have focused on answering it. And some of these earlier attempts one of the things that they struggled with was really high false positive rates. The issue being that if you have a lot of sequences and you're trying to look for a motif that is overrepresented in sequences, you almost invariably would find something. So the robustness of the tools early on was something that was an issue, and over time it has got a lot better as we have been better at solving these kinds of problems and better exploring and explaining and even defining the background population of sequences. Now, when I was a graduate student, I met two postdocs in Saita-Wazir's lab at Princeton who were coming at this from a different perspective, which I really enjoyed, and I learned so much from them, and it has really shaped how I look at the problem broadly. And these two postdocs were Noam Slonim and Olivia Lamento, who have since gone to do many other things. But the problem that they were working on when I met them was actually a solution to this motif discovery problem that was kind of like more universal approach and also quite sensitive and with very low false positive rates. That was, that was really the the goal here. So what was the solution? The solution was really to kind of bring information theory into the fold and use mutual information to look at non-random associations between the presence or absence of motifs, as I will talk about them and really explain what I mean by motif. But, you know, broadly, the presence or absence profile of a given element and its behavior in the data. You know, mutual information is a very simple construct. It can be calculated quite rapidly, which back then was extremely important and still extremely important today. So you can actually go through quite a few elements and do really rigorous search of the space, of the sequence space, to be able to quickly calculate different variations of this and do search and optimization. So that's actually one of the key benefits of using mutual information. Another one is that it is defined both for continuous and discrete data, so you can actually have the same underlying approach to answering whether you're looking at a continuous variable, for example log fold changes, or a discrete measurement like, for example, cluster memberships if you have a larger dataset that you have clustered. So, and I think another aspect of mutual information that's important, we'll talk about a little bit, is also the fact that, you know, you can fit things like conditional information into the mix as well that makes this a lot more kind of like a vibrant path to solving other similar problems as well. Now, to kind of like point it out, here is a few examples. For example, we are looking at 3 clusters of genes, and this specific motif, as you can see, is overrepresented in this cluster. Or similarly, if you're looking at log fold changes, you see there is an enrichment of sequences at the top, so meaning that this motif is associated with higher expression in this data. But what is important here is that Mutual information can also capture nonlinear kind of patterns. So for example, if you have a dual regulator that is simultaneously up- or down-regulating genes, mutual information will be significant for them, whereas with more traditional measures, that might not be the case. So by and large, information theory is a very useful tool for kind of answering this problem, and we have over the years really So we kind of stuck with this as a concept to play with in this space. Now, the way that FHIR— so I should also mention that the tool that does this, that was developed by Noam and Olivier, is called FHIR. And the way that it actually solves the problem is starting to search the Cambridge space. You as the user can define it, but I think they basically started with 7-mers. And what you do is you enumerate all possible 7-mers and you use mutual information to score them. And then you basically sort them very quickly from, you know, highest information to lowest. And then you can do more rigorous optimizations to go from k-mers to— and this is done kind of like in a greedy search— to go from a k-mer representation to a more degenerate motif representation that captures multiple significant k-mers. Now, of course, I should mention that the motifs here are represented kind of as regular expressions, and I'll get back to this as well. So one of the questions that I was working on as a graduate student was actually, again, like early days of applying these kinds of tools and these ways of thinking to the explosion of data that was emerging, gene expression data that was emerging in the cancer field. And we found quite a bit of information in the promoter sequence and what is likely to be transcriptional regulation, but also quite a bit of information, regulatory information, in 3' UTRs, where are known to be sites of post-transcriptional regulation. So I started to pay more and more attention to kind of like that aspect of gene expression control, but One thing that is important to understand is that when we go from transcriptional control to post-transcriptional control, we take some understanding and some outlook on how this happens from one field to another that has really helped in some instances and not so much in others. You know, I was not the only one to do this transition. There have been quite a few scientists in the field who start— who started really their work in gene expression control on transcriptional side because that was really understood to be the main path of kind of regulating gene expression, and more slowly gravitated towards post-transcriptional regulation. And again, we took some lessons learned from that field, including motif discovery as an example, but also, for example, really thinking about RNA-binding proteins being equivalent to transcription factors. And again, this helped jumpstart a lot of early successes in the field, but also because these, these are not like perfect counterparts of each other, as there are some aspects that, that we really needed to step back and rethink post-transcriptional control really on as its own problem. With its own really idiosyncrasies and kind of molecular mechanisms. And one of the key aspects that we needed to capture is really the impact of RNA structure on regulatory interactions. And even at the time, there were many instances of regulatory RNA structures. For example, iron response elements are a good example, or IRIs is another example. There are RNAs that Basically, they're structured and the function depends on it. It's really scaffolding RNAs. There are RNAs like tRNAs that are very well structured. Primary microRNAs are very structured, and that structure impacts their function and how they're processed. Basically, there's this wealth of— even then, there was this wealth of information around the role of RNA secondary structure, or ultimately tertiary structure, on RNA-level post-transcriptional regulation. But for the most part, the way that we were kind of like looking for cis-regulatory elements in that space was really ignoring the role of secondary structure, and that's certainly something that we needed to account for. One way that we can look at how important accounting for RNA structure is can be answered by kind of like a a simple kind of a thought experiment. So one of the, one of the ways that we know about our RNA-protein interactions is through a method either in vivo or in vitro that is very similar to, you know, how transcription factor to DNA sequences binding to these DNA sequences are measured. So these are basically affinity type experiments, and what you can do is basically ask if I have a data set of binding patterns of different RNA-binding proteins, and this is, for example, publicly available data for some of these RNA-binding proteins. If I basically focus on tools that find sequence elements, basically linear cis-regulatory elements, how well do I— would I do in capturing the binding patterns? And, you know, Here I'm showing you two ways of answering these questions for each of these RNA-binding proteins. And as you can see, for some RNA-binding proteins, we do extremely well. It's very clear that they are sequence-driven. For some others, not so much, despite the fact that they have RNA-binding domains. For example, hnRP2B1, it has multiple RRM domains that have a sequence-specific binding. In reality, there is a context to this interaction. And actually, one of the earliest examples that we found that secondary structure actually plays a role in, in the direct RNA-protein interactions was actually in the context of A2B1, where we found a lot of these motifs are embedded in a structured form that A2B1 binds to, and it has a regulatory role in that context. So now, uh, now that I've kind of explained why we want to solve this problem the way that we want to solve it, um, I want to kind of again go take a step back and do comparisons with the linear element, um, discovery side, because that is something that folks have thought about, and maybe that's a good starting point. Um, on the linear side, I mentioned there are different ways of looking thinking about how can you represent a motif and how much information that kind of representation holds. There is no perfect representation. Any representation is a summarization of all binding patterns, so you'll lose some information. So the question is, for your use case, what is the level of abstraction that you're happy with? For example, FIRE kind of has 3 representations. It starts with kmers or seed sequences. Sequences, basically. So these are the simplest forms of quote-unquote motifs that you can work with. Then extends these to regular expressions where you basically allow for multiple options at different positions along the sequence element. And finally, what folks really define as linear motifs are position rate matrices, which offer more resolution in terms of binding preferences. Now, the difference is that, you know, basically position weight matrices is an infinite space. You need some sort of heuristic to get there. Regular expression and K-mers can be actually traversed computationally, so you can cover for the most part that space, and certainly K-mers and less so the regular expression. Now, here we're starting with this understanding that these solutions really work in the linear space, and we can actually find motifs that explain either DNA or RNA behavior. Now, the question becomes, how can we extend this to structural elements? What is the right way of really thinking about representing structural elements? You know, we all understand that RNA structure is a thing. It folds onto itself, and that can play a role, but Whether these foldings are then part of a shared grammar across many transcripts is basically another question, right? So, and how can we actually find them? So one way that everyone knows how we represent RNA structure is this kind of representation, when every sequence is basically accompanied by a representation of this structure and using parentheses or brackets, the idea being that this is a central loop. For example, this is a loop, and then you have stems. And these parentheses basically tell you which nucleotide is interacting with which nucleotide. And together, these two representations allow you to regenerate what kind of secondary structure you're actually looking at. And this is kind of the output of a lot of RNA folding algorithms out there that you put in a sequence and you get this structure. Now, uh, what is important is that, uh, in here, as I mentioned, in the RNA folding field, you're trying to kind of like minimize energy and kind of solve a different problem, biophysical problem, in terms of how do you go from here to here. But in the motif discovery field, we're trying to take many representations of sequences and ask ask whether they have matches to potential secondary structures or not. So it's not necessarily the case that the structure we are looking for is the most energetically favorable conformation of that RNA, if you will, but that doesn't really preclude it from being a functional element. So in our approaches, we kind of stayed away from RNA folding, but as we'll discuss today, that's— you can certainly always bring it in and kind of like add it as a component to answering this problem. So let's talk a little bit about how can we actually represent secondary structures in a data structure that allows us to really work with quickly in the in silico, basically. So one of the key kind of like concepts here that really helped us and the field as a whole to really answer this question is basically context-free grammars. You know, I don't really have a lot of time to go into this. This is borrowed from linguistics. I'm just going to basically present the two versions of this that are useful for today. You know, something like a regular language, you have elements that are constant, and then you can basically link these constant elements together to form longer and longer phrases. So for example, in the regular language, you can create a K-mer or a regular expression even. Basically, the idea is that you start from nothing and you keep adding, you know, C-A-G-U, and then you get a K-mer. So each one is a rule that works on the previous rule and creates a new phrase, and stringing them all together gives you a K-MER. Now, with RNA structure, we want something that can capture not just the sequence but also the information around about what goes with what in the secondary structure space. And the system needs a little bit of memory. So for example, you need— even you see this A, you need to know that it goes with this U, and these two are linked. So if I change this, I need to change that for the structure to remain constant. And so this little bit of memory is actually captured in context-free grammars, where you can basically start with what you had in regular language, where, you know, for example, here we are recreating this loop here by just adding AG. It's kind of a K-mer. But then after that, we can actually add elements to both sides. That captures the fact that these two actually go together. So, you know, U and A go together, and U go together, and G and C go together. And again, running through all of this will recreate this representation for you. So basically, context-free grammars are a data structure, provided data structure, that allows us to very quickly represent kind of like secondary structures and sequences together in one place and really enable kind of like motif discovery in the structural space, if you will.

Honey, quick question. So what was that like, a linear versus polynomial that you had in the previous slide?

Basically, it's basically growth in terms— so basically the length basically grows linearly as you go through, but here it can go polynomial.

Thank you.

So one of the first tools to actually use this— I mean, there's using context-free grammars to represent RNA structure actually goes back to '90s, I believe. One of the first tools to actually take advantage of it to do motif discovery in some sense in the structural space was this paper It really had a profound effect on myself when I was a graduate student and I read this. So I will certainly encourage you to take a look at it if you have time. So it's basically a tool called RNN-Promo where they showed that you can actually take advantage of context-free grammars, and specifically the variation of them that are stochastic context-free grammars, to actually learn what are the shared secondary structures in kind of sequences that, you know, are again like shared regulatory elements. This is very similar to beginnings of kind of the linear sequence space as well, where we were starting with a few genes that we knew that are behaving similarly in certain conditions, and we knew that they need to have similar sequence elements, and that's really how the search started. Here we are doing a kind of tracking the same progression for the RNA structure space as well. So this is really where we started, being able to find these structure elements using and representing them basically using stochastic context-free grammars. So I said all of that to say that basically in the structural element space, we can use context-free grammars to actually represent, you know, any variation of the linear space, like, you know, both K-mer, regular expressions, or even position-weight matrices. In the kind of like a visualization space, you have kind of the sequence-structure duo, which I talk about. You have these you know, secondary structure representations. You can use, again, degenerate sites where you can have multiple options. This is basically equivalent of regular expression as a schematic. And ultimately, you can use context-free grammars to basically represent all of this in silico, ultimately. Okay, so how do we take advantage of this? So the way that I want to think about— I want you to really think about this— is that, you know, context-free grammars can basically play the role of k-mers in the original problem, motif discovery problem, that we were trying to solve. So when k-mers were basically the seed sequences that we used for linear motif discovery, basically context-free grammars can also serve a similar role. The difference, the major difference, is that the search space is substantially larger. So we cannot enumerate all possible seed context-free grammars, if you will. We need to somehow define them, what we mean by seed structure elements. But in some sense, we are starting with some sampling of seed structure elements in the massive space of structure elements. And again, we try to put constraints as as you would imagine, on what this search space would look like. But the idea is that each of these elements will capture some sequence or some structure restrictions. And again, initially, you can start with actual k-mer representations embedded in this sequence. So for example, this is really a G-A-G sequence embedded in this structure. The idea being that if you know these two and you know the structure, you know these two, right? That's kind of how information is shared. And so if you have this, and of course you can further generalize this to kind of like similar to regular expressions where you can have multiple options in each of these positions, but this can be basically our seed representation. And if you can start with a large sampling of these along this space, again, Your computational power dictates how large this space can be, and based on that, how complex these seeds can be. Initially, you know, you can start with like a certain limit on the length of the stem or the loop or how much information is in the sequence content that we'll actually show examples of today. So once we have generated representations of these Cs, again, very similar to K-mers, we can use, again, context-free grammars to scan sequences for matches very quickly and do what we did for K-mers. If we see a match— so for example, here I'm showing you this schematic where we implemented in our tool, which we call TEASER, for structural discovery. It's basically an extension of FIRE that was developed by Noam and Olivier to the RNA structural space. And this was work that I did in collaboration with a close friend of mine, Mohammad Najaf-Awadi, who has a computational lab at McGill, who I encourage you to check out. He's really doing amazing work these days as well. But we worked on this as graduate students together, and we basically developed TEASER to be similar to FHIR but be extended to the RNA structural space as well. So the idea is that once you have definitions of these seeds, you can enumerate them and ask whether they're present or absent in a given sequence. And then you have some measurements similar to gene expression, as we talked for FIRE. Here we are looking at RNA stability, for example. Since RNA stability is something that's post-transcriptional regulated, it's a good example to have. So for example, this element is very clearly associated with destabilization in RNA. And again, you can capture this using mutual information. In practice, the search algorithm is very similar to what happens in FHIR as well. Again, you start with these structural seeds, you use mutual information to capture their association with the data, and then you pick each of them and then you do a more local greedy search around that seed for better, more complex representation. You can go to more complex representations from basically a K-mer representation to more regular expression representation. And along that path, you can also make changes to the structure as well, which again, I think Matt today will give you some examples of. So once you have a more archetypal kind of representative of an element, kind of like the best k-mer, if you will, of a structure element, then you can basically create a bag where, you know, the way that you can think about a regular expression representation of a linear motif is basically a bag of k-mers. So in the same way, you can put together a bag of context-free grammars that look similar to basically create your representation of a structural element in RNA. And so that's basically how our first attempt to solve this problem actually took shape.

So I got a quick question.

Sure.

Can you give us an intuition about, like, how often do you actually see matches? Let me try to rephrase that. So it looks to me that on the one hand, if you choose a context, like a specific, like, you know, context-free grammar, that, you know, the odds of being able to, you know, match it in any sequence is pretty small, right? But on the other hand, you have a lot of, you know, different context-free grammars and you have a lot of different sequences. So there's, you know, that there's an explosion of options in that sense too. So, um, so I'm kind of like not, you know, I don't have the intuition, you know, what's going to be— for example, what's the— what does it mean for your false positive rate?

Good, very good question. So actually, uh, I think Matt will, will talk about this. Uh, it really goes back to which space you start your sampling from. How you define this space is actually kind of problem-specific because you want to capture seeds that are not present in everything and not absent in everything. You can define how much information you want it to have and how many instances you want it to have in that dataset on average that you're starting with. So sequences— so for example, if you have— if you're searching in very long sequences, then these need to be a lot more specific to begin with. But if you're looking for, you know, in shorter sequences, you can make them more and more degenerate from the get-go. So it really is problem-dependent in some sense.

Thank you.

Unfortunately. So I'm going to spend the final few minutes to talk about where we applied this early on, and then I'll, you know, for a second talk, we'll hand off to Matt and Mehran to talk about where we are with kind of this, with this problem today. So as I mentioned, I'm interested in metastasis as a scientist, and one of our go-to models for answering questions about metastasis are these really nice in vivo selection models where we basically start with established cell lines which we call the parental line that are often poorly metastatic. One of our favorite models is actually the MDAMB231 model. That was generated by Joan Massague's group at Memorial Sloan Kettering. Basically, they started with this cell line and then they selected this line in mice for higher metastatic potential and they derived, for example, the MBLM2 model, which is substantially more metastatic to the lung. LM2 stands for lung metastatic round 2. So, basically, 2 rounds of selection in the lung gives you the LM2 model. And And, you know, both the Massague lab and other labs, including ours, have over the years generated a lot of these matched models of metastasis. And because they are really isogenic but phenotypically very different, it is a very good model, matched model, to really explore molecular drivers of metastasis. And what is really important is that even in this, like, short selection going from the parental to the LM2 models, Again, you see broad changes in the gene expression landscape of the cell. Again, like, even a short selection like this really brings about all the regulatory modulations that you need to see. So, of course, there is, like, a ton to learn from these models in terms of gene expression control. But at the same time, since they have more or less that, you know, they're isogenic, it's also very, you know, nice control background in that sense. So, uh, one of the, one of the first, um, kind of like measurements that I did as a postdoc, um, was to use a metabolic labeling of RNA to measure RNA decay rate in these models. I actually developed this, this method, um, in, in grad school. There are quite a few groups that around the same time, uh, did similar measurements, uh, that with the goal of trying to measure RNA decay rate transcriptome-wide. And the idea is that you feed the cells like thiuridine, and, you know, which gets more or less seamlessly incorporated into RNA. And, uh, you know, this is basically the pulse period, and then you remove thiuridine from the media, and then you can basically, through a chase period, measure how much thiuridine you see in a given RNA, and that tells you kind of the stability of that RNA. So that's basically the idea. So I did this experiment with, uh, prolylhydraminostatic cells in biological quadruplets. So we actually had really nice data comparing poorly and highly metastatic cells. And this is the results of that data. So, again, quite a few genes that are up- or downregulated as you go from a poorly metastatic state to a highly metastatic state. And here I'm showing you two major clusters, as you, you know, it's comparison, so there's really two main clusters you're interested in. Those are stabilized and highly metastatic cells, meaning that they have like lower decay rate, and those are destabilizing hyaluronic static cells. We had some intuitions about what's going on here. My, my postdoc advisor, Soheil Tavassi, actually had done a lot of work in this area. He was one of the pioneers of looking at microRNAs as regulatory factors in cancer, and we had some intuitions about like what, what factors could be responsible for this, but not not much known was about this, this set. And this is what we tried to answer, not necessarily using structural element discovery, but we looked for both for linear and structural elements using FIRE and TEASER. And out of that came this, what we call basically RNA structural element, structural RNA stability element. Sorry, the S stands for structured. Here I'm showing you this element as a generic stem loop. This is just to show that it's a fairly GC-rich stem loop element. And what TEASER actually shows is one representation of this along with the mutual information and z-score, but also gives you— generates a heat map that tells you where things are actually enriched or where things are depleted. So, for example, if you look at transcripts that are destabilized with respect to those transcripts that have these SRSC elements, we see a very significant overlap. You can actually calculate the p-value for this using hypergeometric distribution and, you know, define a minus log of this p-value as enrichment score and use it to draw a heat map. On the other side, in the background, for other transcripts, SRSC overlaps are actually depleted, and that gives you an enrichment score of -8. And, you know, you add that to your heat map. So basically you get this like heat map representation of where these elements are along the input data as well, in addition to the fact that it's just a non-random pattern of presence or absence. And to give you some instances of these elements that we find, again, very GC-rich sound loop structures and that are overrepresented in these destabilized elements. What I showed you in the previous slide was basically in the kind of like discrete space. I'd group them into destabilized and not, whereas in reality you can actually run this in the continuous space going from stable to unstable with substantially more bins. So each of these bins have about 350 transcripts in them, and again, as you see, there's a very nice significant enrichment of these SRST elements on the left, and they're depleted on the right. So this is exactly the kind of pattern that we would expect to see for functional elements. And of course, you can, you know, if an element is associated with RNA stability changes, you would expect it to also be associated with gene expression changes, and that's exactly what we saw. That this element is enriched among the genes that are downregulated. We can look at it in this heat map context, but you can also look at it— basically look at the density plots. Again, this is the background, and this RSC regulon, defined as a set of transcripts that carry these elements, they show a significant shift to the left. So this, this kind of answers the first question that we always want to answer when we go from transcriptomic measurements to some understanding of what the underlying cis-regulatory elements are. The second part is also important. So what are the transactors that are actually talking to these elements really forms the regulatory mechanism of gene expression control. Over the years, we have developed a couple of ways of answering this, either computational or experimental. Experimentally, you can basically do variations of bait-IP mass spectrometry, or more recently genetic screens, to be able to find factors that are actually interacting with this, with these elements that you have found. Or you can actually do this computationally, or at least attempt to do this computationally, and this is based on, again, co-expression analysis. The idea being that if you have a regulator that sits upstream of a regulon and your regulon has a given expression pattern or RNA stability pattern in this, this example, you would, you would expect the regulator to be co-expressed if it's an activator and anti-correlated if, if it's a repressor. This is a very simple paradigm, it doesn't always hold true, but sometimes does. Actually, in my experience in the RNA binding protein field, it's a, it's a pretty good proxy and starting point. So this is where we start usually to to be able to find trans factors, but not always successful. So for example, if this is your target regulon, you'll find RNA-binding proteins that are correlated or anti-correlated in terms of their patterns in existing datasets. So we did this analysis for the SRIC regulon that I showed you was associated with RNA stability, and we found 3 potential RNA-binding proteins that, that could be responsible for for the function of this element. And one of them we could actually validate using knockdown profiling. So we found this protein called TORBP2, which when you knock down, you actually see an upregulation in this regulon. Again, you can also look at the density and the shift to the right in this case. And this shift is actually due to higher stability of this regulon, as you would expect, because we are looking at stability element. And if we overexpress TORBP2, we see the opposite effect. So that's very clear, kind of like one-to-one correspondence between TORBP2 and activity of this element. And finally, to really prove that there is this direct interaction, you can take advantage of a method like CLIP that tells you where the RNA-binding proteins are sitting transcriptome-wide. So for those of you who don't know CLIP, it's basically equivalent of ChIP-seq for RNA-binding proteins, it actually relies on UV cross-linking, so it's a lot more precise and gives you like nucleotide resolution of RNA-binding protein RNA interactions transcriptome-wide. But you can basically look at where TARGP2 sits and ask whether they actually form these structural elements that we were looking for, and the answer is by and large yes. So TARGP2 clearly binds elements that that are similar to our definition of SRSC elements. And this wasn't actually too surprising in the sense that TARBP2 is actually a double-stranded RNA binding protein. But this role that we found for it was actually a non-canonical role for this RNA binding protein. In its canonical form, it actually interacts with Dyson and Argonaut and is a component of a RISC complex machinery. But we also showed that it directly interacts with RNA as well and has a regulatory function. And we had a follow-up paper a few years ago that really looked more closely at this problem. But what I wanted to say is that, like, all of, all of this pathway of metastasis that is modulated by TARBP2 that we found really for the first time in the context of breast cancer metastasis was really initiated by teaser, really, how we went from transcriptome-wide measurements of RNA instability to kind of like a new understanding of how gene expression is regulated in cancer cells, but also just broadly. So this is basically where things were around the time that I started my lab, and over the years there have been a few changes in the field that has also modified how we are approaching the problem. And these are basically what are going to form the foundations of the next talks that you will hear about. One thing to point out is basically in recent years, thanks to advances to RNA structure probing, we actually now have tools that tell us what RNA looks like in vivo, even certainly in vitro, but even in vivo. And there are a number of transcriptomic tools that answer these kinds of questions. For example, DMS MapSeq, that was developed by initially a version of it by Sylvain Ruskin in Jonathan Weissman's group, and IC-Shift by Ryan Flynn and Howard Chang's group. They basically— both of these work based on RNA modifications that are a function of how exposed the RNA is. So they basically label different nucleotides in RNA based on whether they are part of a loop or a bulge or part of a stem. So you can use them to kind of get a sense of the structure of the RNA. And that is basically what prompted the project that— in part prompted what the project that Matt did, basically that turning TEASER into piTEASER. That can also incorporate other sources of information as well into not just sequence information as part of the search space and also really just make something that's a lot more, you know, Pythonic and modular than what we had initially. Another area that really the field has gone into, the entire gene expression regulation, is really the advances in deep neural networks and how they are reshaping how we model gene expression control and cis-regulatory elements, you know, starting from, you know, DeepBind that I showed you some results from earlier, but larger and larger models, both in terms of parameter space or even input sequence, are being used to look at gene expression control both at the transcriptional level as well, but also at post-transcriptional level. So what Mehran will talk about today is that how we try to bring context-free grammars to basically deep neural networks by giving the model basically the rules and let it learn the grammar, basically. And he'll talk about that. So this will allow us to basically take what is now the state of the art for defining regulatory elements not necessarily based on motif representations, but based on model dissections. And there are quite a few tools out there that, you know, including the ones that Anshul Kundaje's lab at Stanford have developed, that you can actually use to interrogate these kinds of models. And from these models, you can actually learn what are the sequence preferences that drive the data And we are trying to— we have tried to bring this into the kind of like this RNA structure dimension into these kinds of processes as well. And I think with that, I'm going to stop here. And I acknowledge both Mehran and Matt who are going to talk today, but also everyone else who has been in my lab and also kind of contributed to these projects.

Hany just gave a great overview of how he approaches the problem of finding cis-regulatory elements by using both computational and experimental methods. And I will elaborate on one angle of this that we implemented in the PyTeaser project. As Hany mentioned, it's a development of his previous work that started with FHIR with k-mer discovery, then Teaser with structural seed discovery. And so this is the current state-of-the-art thing we have in the lab.

To—

give me one second— yeah, to specify the problem, our goal here in this case is to identify RNA secondary structure motifs that explain differential measurement of interest. So Hany used mRNA stability a lot as an example, and I will be using splicing patterns, but really the goal is to make it generalizable framework that can take any differential measurement. The challenges associated with this are that the space of possible RNA-seq construction motifs is extremely large, and our understanding of these cis-regulatory rules are pretty vague. So to give you a framework of what our solution is, the program PyTeaser takes any differential measurement between two conditions as an input, complemented by the set of sequences that we used for these measurements. And it outputs a ranked list of RNA secondary structure motifs that could explain this differential signal. You can also optionally provide RNA secondary structure probing data that Hany has mentioned and that I will talk about more later. So you have a pretty good idea of overall how this program would work conceptually, so I will be focusing more on the details of implementation and the decisions we had to make along the way than on the conceptual part. So I'll start with the seed selection. So as you know, the overall space of RNA secondary structures is pretty much infinite. There are any number, any kinds of them, and since we want to actually try some number of them one after another, we have to put some constraints. And the most obvious constraints are size constraints. So by starting with just short seeds of between approximately 8 to 20 nucleotides long, we can reduce the search space to already 10 to the 11th seeds, something like that. And then additionally, we put information content constraints which I'll talk about more in the next slide. And this way we can actually start with a coarse-grained search using approximately 1 million of short informative seeds. So to go to algorithm of how we choose seeds, we basically just generate all possible seeds of a given size and we filter them based on the information content criteria. We have two key criteria implemented. One of them is number of meaningful bases, so the ones that are not Ns, And the second criteria is restrictiveness of the seed, so probability that sequence of the same length matches the seed. So what I mean by that is the seed can be too permissive or too restrictive. If you take a seed that only consists of n bases, it describes any single sequence of the same length, so it's not specific enough and doesn't really tell us anything. On the other side, if we take a long seed seed where every single nucleotide is defined with specific base. This seed can only cover a very small number of sequences, so it's too specific and doesn't cover enough of the search space. So we want to stay in between, so start with seeds that are somewhat informative but not too vague, and this is what this information content criteria is addressing. So with the default parameters we use, we start with approximately 1 million of seeds. And the package provides Python script for generation of such seeds and also a set of precompiled seeds that we use by default in most projects. I put URLs here because I believe the presentation will be shared after the talk. So anyone interested, feel free to visit these links. Once we generate—

A quick sort of question for you. For the people who are just joining us for the main talk, and you might not have been in the primer. Could you just, you know, in 30 seconds say, you know, what do you mean by a seed in the context of your problem?

Yes, I'll go back here. So I will use two terms. I'll use the word seed and the word motif. So the word motif is some kind of definition of cis-regulatory element. The word seed is we literally start with just a small secondary structure which consists of of a stem and a loop, or maybe two stems and two loops in between them. So these we use as kind of extension of k-mers to identify these structured motifs. And then in the end of the program, I'll show that we group similar seeds together and we call a group of them a motif. Thank you. Yeah. Okay. So once we've chosen a set of seeds, we go to the step where we want to search, to do search, identify matches of each individual seed in all the sequences. This step is very slow because we have a lot of seeds and a lot of sequences. So in order to do it in a computationally feasible time, we had to come up with some high-performance solution. In this slide, I'm just showing you the tricks we use to achieve that. So the first ingredient is context-free grammars that Hany has mentioned. Context-free grammars is a great way to encode, to describe secondary structure for computer, and also to make it possible to search for RNA secondary structure quickly in a sequence. The second ingredient we use is Numba. It's a compiler for Python functions, so instead of running the whole code as default Python interpretive code, we precompile certain function classes, and this can provide us with a large speedup. So on the right, this is not relevant to the specific problem, but I'm just illustrating that Numba, this compiler for Python, can increase speed of computation by more than 2 orders of magnitude in certain cases. And so to specifically implement such seeds, to implement the seeds, we use two core classes. So Numba as Python compiler is great for fast checking but doesn't have wide functionality. So we use one class that only can do one thing. It can search for a seed in a sequence, and it does it very, very fast. And we use a second class, which is this layer of interface between the user and the rest of our program to this Numba class, uh, that is written in regular Python. This class can be used for input/output printing, visualizing, and doing all kinds of interface operation with seeds. So here I'm talking about all these matching of seed and sequences, and I just want to show an example here. So here there's an example seed. It has structure shown here and the sequence shown here. I want to emphasize that sequence does not only have 4 specific nucleotides. It also has some degenerate nucleotides shown by N, which provides larger number of sequences the seed can cover. And here there's an example of a match. So this sequence can potentially form a secondary structure represented by the seed while also matching the sequence. So now we have a number of seeds, we have some initial set of sequences, and we want to do the search for seed matches in these sequences. And if we show each seed here as a row, and then within each row one cell is a sequence, we can show sequences that match the seeds by black squares. And this makes these match profiles that are just binary vectors for each individual seed and whole set of sequences we use. When I say set of sequences, it can be transcriptome, it can be set of exons, anything really. I can say transcriptome for simplicity. So again, at this step, we are trying to find sequences that match the seed of interest, and we do it for each individual seed among the starting set. So just to emphasize again, this step is the slowest part because we have a lot of sequences, a lot of seeds, and that makes the program run slow. The good news is that pretty often a research lab studies mostly one organism, and in this case seed profiles, these matching binary vectors, can be reused between different projects. So we really wanted to add this functionality where we can run the program once for a certain set of sequences and then reuse use such match profiles in the future for other projects. The one problem though is that since there are a lot of seeds and sequences, if we just write these profiles in any conventional way such as tables or anything else, they take a lot of space on the computer. It's not really feasible. So we had to come up with a method to encode the seed profiles in a heavily compressed way. And here on this slide, I am showing how we implemented In the middle here, there's a match profile kind of like the ones I've been showing before. For any single sequence, the only information we really need to know is if it does or doesn't match the given seed. This information is just a binary variable, yes or no, and therefore it can be encoded not even by a byte but by just one bit, and therefore we can encode this vector for n sequences in n/8 bytes. Now, since we have a lot of these profiles, some files might get corrupted and we want the machine to know that, so we add two additional things. First, we add 4 bytes with profile length, so how many bytes does the program need to read, and in the end, we write MD5 checksum and additional 16 bytes to make sure the profile was read correctly. So such thing— so yeah, so we implement this with a Python struct module, and then the size of such individual compressed profile depends on the number of sequences used. But for example, for human genome, which is about 20,000 genes, this compressed profile takes about 2.5 kilobytes. So we can save profiles for 1,000 seeds in only 2.5 megabytes, which at this point makes it feasible. So this is again mostly for the purpose of running this program once for any given organism or set of sequences and then reusing these match profiles later. So now I want to talk about the NuEvo probing. So I've been mentioning this searching for seed match. So when we say a seed matches a sequence, we say that sequence can potentially form secondary structure that is described by the seed. We don't actually know if it does form a structure, but we hypothesize it can. We can actually check if it does form this structure in vitro, in vivo, or in silico. And as Hani mentioned, within the last 5 years, these methods for in vivo for secondary structure probing have been developed. The two most popular methods are DMS-Seq and SHAPE. They are pretty similar conceptually. In this slide, I'm showing the scheme for DMS-Seq from the original paper. You use this highly toxic chemical, DMS, that preferentially binds— well, preferentially interacts with the flexible or rather unpaired nucleotides in the RNA. And then you take this RNA and you proceed with the library preparation. And depending on which reverse transcription enzyme you use, DMS doesn't really help this reverse transcription enzyme do its job, and it causes either truncations or mutations of the reverse transcription process. So GNA will be either mutated or truncated. We use these mutated or truncated reads for further library preparation and sequencing, and we can count them in the end in the sequencing data. They're shown here by these red dots, and basically the data in the end looks like this profile where for each nucleotide you have a number how often this nucleotide is mutated, and that shows you approximately how often this nucleotide is unpaired within the population of RNA species. The second popular method that's called IC-Shape is very similar conceptually. It uses a different molecule to do this probing. This molecule is less toxic. In the end, the data you get is conceptually pretty similar. It's also a proxy for how open each nucleotide is, and it comes as this profile, one number per nucleotide. So how do we apply this data to our motif search? So we found a match of a certain seed, and we hypothesize that the sequence can form secondary structure we described. Now we pretty much take the DMS or shape DMSC data and overlay it on top of this sequence. And here I'm showing it in these colored nucleotide circles. The darker a nucleotide is, the more mutations we have in this DMSC data, the more open or unpaired this nucleotide is in the population. So here you can see that actually the nucleotides that we predict as unpaired are generally darker than the nucleotides within the stem. So this does not contradict our hypotheses that this sequence forms the structure shown here. So we call such comparison a pass. So when we go back to these whole seed match profiles, we can use DMS or shape reactivity profiles to filter out the matches that are possible just based on the sequence but unrealistic if we look at in vivo secondary structure data. So where we see the conflicts of the predictions and the in vivo data. This does not have to be restricted to just in vivo profiles. We can use in vitro profiles. We can even use the in silico secondary structure folding, which is a less favorable option but currently far more accessible than in vivo protein data. So at this step, we filter out the matches that are unrealistic. And now with that, we go to a step where we discretize the expression profiles and calculate mutual information. Han has shown it with this beautiful scheme on the right. We start with some kind of input measurement that can be expression stability, anything really. I visualize here by a volcano pot, would discretize such measurement by bins shown here and end up with a discrete measurement vector. And this is really the point where we can estimate how good a certain seed is for explaining our differential signal. So we calculate mutual information mentioned before between the seed seed profile and this discrete vector of measurement of interest, and this mutual information value is our estimate of how good the seed is in explaining the signal. So, we can do this fairly quickly for all the seeds we started with, and then we can rank the list of seeds by the mutual information value, assuming that the higher the seed is in the list, the more likely it is to be explanatory. Obviously, we don't want to proceed taking all the seeds within the list. We want to— oops— we want to take some limited set of significant seeds, so we use a binary search algorithm. So for any single mutual information value, we can calculate the p-value associated with it, by doing permutation tests, calculating mutual information many times for shuffled vectors, and based on these permutation tests, we can search, do binary search for a good threshold with a certain value of FDR and only take the top number of seeds with the highest p-values, I mean, lowest p-values. So at this point, we have a small number of significant seeds. Once we have a significant seed, we want to optimize it, and the logic behind it is that in the very beginning I mentioned we do a coarse-grained search of the whole RNA secondary structure space. We can't possibly do a fine-grained search, so whatever good seeds we find are probably suboptimal. There's probably some local maxima close to the seed, so we're trying to find this local maxima by iteratively improving the seed. So we take any single nucleotide within the seed and consecutively mutate it and try to find modified versions of the seed that have higher mutual information value. We don't have to only change the nucleotides within the seed. We can also add nucleotides to the end of the seed, so prolonging it and adding information there, basically, and also do it in an iterative way, trying to find seeds with higher mutual information value. So, every single seed within this short significant list is being optimized this way, and once we optimize it, we often observe that these seeds that are highly optimized are actually— many of them are similar to each other. They have similar structure, similar sequence, and similar profile, and they probably represent the same rules in slightly different ways. So it makes sense to cluster these seeds together, and so here is where the word motif comes. So here we call a cluster of seeds that are similar to each other a motif. So such motif is the final output of the program. For any single set of differential measurements, we identify a motif that could explain such difference in the signal. So now I want to talk about one application we use it for. We use TEASER to study splicing in metastatic cancer. As Hani mentioned, metastatic cancer is one of the key model systems in the lab. So we analyzed cassette exons that are differentially spliced between highly and poorly metastatic cells, and we used PyTeaser to identify motifs that would explain such difference. The top motif we identified is this motif shown here. It's very CG-rich element. We termed it SSE as structural splicing enhancer. Here there's a diagram that shows that this motif is highly enriched on one side of the spectrum if we rank the exons by their relative splicing score. We wanted to know what transfactor could bind this motif, so we performed an RNA pulldown. We used an RNA that contains several copies of this motif, and we followed this by mass spectrometry. And the top hit in the mass spec samples was SNURPA1, Snorkel1, it's a part of spidosome, but it's not known to bind anything like the element we found. So we wanted to know if Snorkel1 actually has preference for motifs like SSE. So to do that, we performed RNA binding-seq in vitro experiment. It's an experiment where you basically mix a bunch of RNAs together and then sequence, trying to see which RNAs are bound by the factor preferentially. So on the right here there's a diagram where what I wanted to show you is that we tried to perturb the motif in different ways and see how it changes SNRPA1 binding to the motif. So if we— so there's a wild type on the left here. If we modify the sequence in a way that preserves the structure, for the most part the SNRPA1 binding preference stays the same. If we disrupt the structure, the binding drops significantly. So you might think structure matters and sequence doesn't. It's not quite the case because if we also perturb the loop sequence, we get even worse drop in binding affinity. So the outcome here is that basically SNRPA1 has strong preference for motifs like the one we found, and both sequence and structure matter in this case. We proceeded to validate the STERPO1 binds structures that look like SSE in vivo by doing cross-linking and immunoprecipitation experiment, which is pretty much a gold standard in this field. And further, we validated that STERPO1 is driving cancer metastasis by doing mouse experiments where we knocked down SNRPA1 and measured the cancer progression and observed significant change with SNRPA1 knockdown. So the message here is that we used PyTeaser to identify non-canonical splicing regulatory code in metastatic cancers. This pathway is orchestrated by SNRPA1 binding structured motifs and changing the splicing patterns. And that's about it for the application. So the conclusions for my talk are that PyTeaser identifies RNA secondary structure motifs. It employs context-free grammars and mutual information that Hany emphasized in the introduction. And we used PyTeaser to identify STERCO1 as a non-canonical driver of metastasis. So with that, I would like to thank especially Hany me, but also everyone in our lab for creating a very supportive environment. And I'm happy to take any questions.

I will be talking about using this concept of context-free grammar for predicting RNA-binding proteins by using neural networks. So basically, there has been quite a few methods that earlier Hani was pointing at them, and the contribution that we want to make here is to use convolutional networks in a very specific way to capture the context-free grammar and then provide it for the later layers of the neural network to learn from them. So we have three aims that I'm summarizing the project with. Initially, I will talk about how we can use this specific type of convolutional layers to make sure they can capture both sequence and structure.. And then once I've shown that, I will show how this network is capable of learning RNA-binding protein binding sites as obtained in vivo from CLIP experiments. And at the end, we will try to dig into these learned networks and see if we can interpret them further. So the name from the project comes from this one version of this Greek mythology about Pythia, where she was the priestess at the Oracle of Delphi. And when she would go to this water spring to drink from it, sometimes there was this python that would show up and come to this water spring and start to communicate with Pythia. So when this communication was happening, Pythia would go to this ecstatic state, and she would make these unclear sounds. And when the priest was passing by, and when the priest was hearing these sounds that at the beginning they seemed meaningless sounds, actually he was able to interpret them and find out that they were actually meaning the answer to life, the universe, and everything. And the analogy that I want to make here is that the python coming to the water spring is basically an RNA RNA binding protein. And PyTHEA is our neural network which is seeing this interaction and is so amazed by it that it goes to this ecstatic state and learns all these millions of parameters. And of course, it will be our job as the priest to really see what these millions of parameters mean by trying to interpret them. Basically, I will go to aim 1, how we can use convolutional layers in a specific way to make sure they learn a structure. Basically, let's start with this particular sequence that I have made sure it represents a short stem loop. All the bold nucleotides that you see perfectly base pair with each other, and we have this 3-nucleotide AAA in the middle. And if you want to design a convolutional filter that can capture a C and G that occur with 3 nucleotides from each other, that would be very similar to a convolutional layer that has 4 input channels representing the 4 nucleotides, 1 output channel. It has a kernel size of 2, meaning that it sees 2 nucleotides at a time, but also it has dilation of 3 meaning that it will look into 2 nucleotides while keeping 3 of them. So with this concept, we can define the weights. So when C and G occur, a specific number occurs, and the other cases, other values. And we want to just detect this one C and G base pairing in this structure. So if we start by shifting this filter across the sequence and calculate the cross-correlation, it will produce different values at different nucleotides. But when it gets to the occurrence of CNG that is specifically 3 nucleotides apart, it will produce its maximum value, which in this case is 2. So this is like one convolutional filter with dilation of 3, and we can do this for all the other dinucleotide base pairing that we are interested in. So for this specific sequence here, from left to right, I'm showing all the possible base pairings that occur to make this perfect stem loop. And you might think that, okay, and now that we know convolutional layers can learn this, why don't we just let the network learn them? And the idea is that basically if we want to have like one convolutional layer, we need to define what is going to be the dilation. And if we can manually define these weights, we can use zeros in a way that would reflect all the possible dilations. So instead of like passing the sequence multiple times through multiple dilation, you can do all of that in one pass. So for this specific sequence, I'm showing you how the output of these manually designed filters will look like if you have a perfect ISLM loop. So at the center of the left sequence where we have a perfect ISLM loop, you see that the highest values align with each other. But if you permute the right side of this sequence, the structure becomes different than the perfect ISLM loop, and you will see the maximum values are at different places, and you are not capturing this. That is specific design. So here, basically, the goal of this manually designed convolutional layer is to provide a context-free grammar. It's not learning all of the structure. It provides these numbers, and then it's up to the next layers of the neural network to map this to the, basically, objective function that we have at the end, whether it's a structure or binding of a specific RNA binding protein. So to show you why we need to design this fixed dilated layer as opposed to passing the input sequence through all these different layers, I did a toy experiment and I designed this dataset that I'm going from the length of sequence that is 50 up to 1,000, and I'm engineering my fixed dilated layer so it's going to represent basically dilations happening from 4 to 50. So I'm not even going to the full range of dilations, just a limited range of them. And you can see that the number of minutes that it takes on a GPU to perform just one epoch, it's quite different when we apply this fixed state of the layers versus when we want the model to learn all these basically weights by itself. And this is just one epoch. You can imagine, in a lot of cases, you would want to perform tens of epochs for your model to converge. So to show that this context-free grammar is actually providing enough information for learning some aspects of RNA structure, I'm using the minimum free energy of RNA sequences as a surrogate of RNA structure. So my task here is to show that that when I provide this input sequence to the network, at the beginning it goes through this fixed dilated convolution that its weights are fixed, it's defined by us based on prior knowledge of base pairing. And then this context-free grammar captured at the fixed dilated layer is passed through the other layers of the network so they can map it to a minimum free energy. So I did this for a few different lengths of sequences from again 50 to 1000, and on the top I'm showing you how the loss was converging as we were performing our training, and also how the Pearson correlation of prediction versus the minimum free energy of folding was looking like. So So for the on the left that we have a length of only 50 base pairs, it's an easier task. But as we get to the right where we have length of 1,000, it gets harder for the model to learn these patterns. And you can see that, for example, at the beginning, in just a few epochs, the model really converges. But on the right where we have like a long sequence, it takes a large number of epochs for the model to converge. But in all cases, the model converges quite well, and at the end you have the same Pearson correlation on the held-out dataset. So, so far this was just focusing on one task, that this is capturing a structure, and what we want to do next is actually take this and bring it to the context of predicting in vivo RNA binding protein sites. So for that—

yeah, good question. So, so the minimum free energy that you used here, is it, is it from some, um, some simulation or is it experimentally measured?

It's, uh, yes, so it, uh, it's actually from, um, one of the, one of the, uh, these, uh, chemodynamic simulations, uh, from VNRNA. It's not experimental data.

Okay, and then, uh, yeah, thank you. And then in the, in some of the fixed dilated convolution layer in the beginning, like, um, how many different sort of seeds, so to speak, are you including? Is it just like one type of structure or like many, many different structures?

So basically, um, here for the fixed dilated layer, for this specific task, I was using, uh, the dilations that range from, uh, 2 base pairs to 96 base pairs, and in the case that the sequence was longer than, of course, 50 base pairs. So I was doing— I was not going above like 96 in most cases. And in addition to like dilate— fixed dilations, I was also including bulges. So meaning that if there's some asymmetry on one side, and that frequently happens on many RNA structures, it can also be captured.

Okay, thank you.

No problem. So basically, so far we have some evidence that at least based on thermodynamic simulations of RNA, of that map, and identify like the best folding for RNA PyTHEA, or at least the version that I was showing you so far, Pythia was able to map sequence to its minimum free energy. And now we want to make one modification and try to predict— try to see how Pythia can predict antibody proteins as determined by CLIP assays. So here we made another modification. So while the sequence is passing through this fixed dilator layer, it's also being passed through this regular 1-dimensional convolutional field, convolutional layer. The purpose for that is that we already know that in addition to a structure, sequence is also important. In some of the previous slides that Matt and Hany were talking about, we know that even some specific sequences in the loop part, although they are not affecting a structure, they're quite important for binding of RNA binding proteins. So we have one layer capturing context-free grammar structure and another layer capturing a sequence, and then all of these go through next layers and are concatenated with each other. And at the end, they try to predict a given RNA binding protein profile in a binary classification task. So here I'm showing you the data from one of these RNA-binding proteins, RBM6, and I'm showing you 4 different methods basically. And I'm showing you the receiver operating characteristic curve where it shows false positive rate on the x-axis and the true positive rate on the y-axis. So first of all, this is basically not the ideal way to show binary classification task is specifically in genomic assays where we have imbalance in the data. But we will start from here, and the baseline is this diagonal dashed line. And if you take the DeepBind model that was trained in 2015 on in vitro data, of course we expect that it's not going to work on in vivo clip data that are obtained in in 2021. So if you want to fairly compare with DeepMind, it's better to retrain the model and see how it looks like. And also, I use this other model GraphProt that is using a support vector machine and tries to include both sequence and structure to learn RNA-binding protein preferences. And again, with GraphProt, either trained on the exact same data and tested it out on the exact same data. And here we see that PyTIA is doing 6% better than a retrained DeepMind. But on the precision-recall curve, which is a more appropriate comparison, especially in genomic cases where we care really about precision And here we see that there's a 21% improvement from PyTHEA to retrained DeepBind. So this is not just one RNA binding protein. If I look into this through all the 43 RNA binding proteins that we tested, we commonly see the same thing. So these plots that I'm showing you here, the blue bar is showing area on the receiver operating characteristic curve, and the red bars show the area under precision recall. Curve. And PyTorch, I'm highlighting it with purple. If there's GraphProt present, it's with green. And retrained deep bind is specifically specified. And if it was a deep bind model from 2015, it's written without the retrained definition, and it has the motif ID that was used.

So, uh, I'm a little bit confused here. Sorry, uh, So, so, so you know, your neural network is only taking the RNA sequence and it's predicting that it's going to bind or not, right? But it's not taking any representation of the protein that is supposed to know whether it's binding to or not.

Is that correct?

That's right. Yeah. So then how does it actually do anything meaningful on held-out data? Because you would naturally think that RNA protein binding is heavily dependent on, you know, the epitopes on the protein and the folded structure of the protein and whatnot, right?

So basically that is our objective function. So based on some— based on experimental data, we train the model in a supervised way, so it sees how a given sequence maps to the binding profile of that specific RNA binding protein. It's not provided as an input, but is the— in a basically supervised training task, it's the output that is being trained on.

Does that explain? Yeah, not actually. So I guess my, you know, let's say you train your model and the first, you know, let's say you, you hold out like RBM6 as you have here, you train your model on everything else here, right?

Um, oh no, no, no, no, the The held-out is not RNA binding protein. So in each RNA binding protein, you're doing a classification task separately. So each of them has their own dataset.

Okay.

So the held-out is not going between RNA binding proteins. It's within the RNA binding protein, you hold out sequences that are—

Sequences.

Thank you.

Sorry. Thanks. Yeah. So basically, PyTorch is generally performing better, particularly when you look at the differences in the red bars, which is the area under precision-recall curve. And one thing that I noticed, as you see, for example, in the case of RBM5 happening here, is that with the default parameters that we used, basically the default DeepMind wasn't able to converge on this state. So Originally, for these models, we used dilations going from 2 to 48, and two of our RNA-binding proteins, DDX5 and RBM42, even with PyTIA, they weren't able to converge. But when we increase the dilation range for these two going from 2 to 96 instead of 2 to 48, then the model was able to converge and provide the same type of accuracy and generalizability, meaning that it's really this dilation concept that is really making an impact here. So, so far our PyTHEA model has this knowledge about binding of each RNA binding protein. So we have one model trained on each RBP, and now we want to obtain some wisdom from this knowledge that PyTHEA has. So here I want to show you, like, how the salience of these convolutionals look like. So what I did here was to perform an in silico deletion experiment, meaning that once I measured the output of the convolution for a given sequence— so here each row is one sequence— and I started by removing 10-mers across this 2056-mer.. And then I measured how the output of the convolution varies from the deleted version versus the reference version. So on the regular 1D convolution, which is capturing sequence, you see that always there's one part of the sequence that is really important and the other parts have a more similar score, meaning there's like this unimodal distribution and the convolution is focusing on a specific part.. And with the fixed dilated layer, which its weights are predefined and fixed, we are seeing this basically bimodal distribution, and it's actually capturing two major ranges of variations when we basically remove 10-mers from the sequence. So I think this provides basically a heterogeneity for the model that it can use this heterogeneity later as its context-free grammar and build upon it to really learn about RNA-binding protein preferences. So to do it the proper way, I used two methods that are developed by Anshul Kundaje's lab back to back. And the first one, DeepLift, is actually using quantifying the contribution from all of the neurons in the network. So instead of just looking at a specific layer and the weights of a specific layer, it's actually quantifying contribution from all neurons. And it's doing that per nucleotide, per sequence. And so we can obtain a set of reference based on the regions that bind a specific RBP and the set of basically a shuffled control and use it with TF Modisco, which is tasked to obtain a robust and non-redundant set of motifs from these importance scores. So I did this for our set of RNA binding proteins, and once I used all of the neurons from Pythia, which is the green bars here, and I also did this once for the 1D once for the regular convolution and once for the fixed dilated convolution. And here on the y-axis, I'm showing you the number of motifs that we discover. Of course, you would expect the motifs that come from all of pythiol together to be more robust, which is the case. But also here we are seeing the extent of heterogeneity that the fixed dilated layer is capturing and providing. So very high number of non-redundant motifs that, and then the later layers based on the non-linearity, would decide to penalize or prioritize for a specific RNA binding protein. So I'm showing you here example from one of these RNA binding proteins, SFPQ, and here are all the motifs that we identified by using PyThL and DeepLift. Here they are clustered based on their similarity according to Continuous Jaccard. Although they're expected to be non-redundant, but still there are clear similarities. But still we can identify motifs from different clusters and compare them and see that there's a general heterogeneity that Pythea is capturing. Some of the parts of these motifs actually reflect what's been discovered before, So, Ray et al., basically they use in vitro assays to identify motifs for the RNA binding proteins. And actually, in many cases that we have the shared RBPs, the motifs show similarities, but also there's some new information that is captured with this in vivo data. But at the end of the day, like in the old days, we would we would discover motifs from in vitro data so we can use them to discover incidences in basically in vivo. But now that we have this model trained for in vivo data, we can actually use it for hypothesis generation in many cases. So one example could be any type of mutation that we want to see its effect in RNA binding protein binding. So for example, if there's a polymorphism that occurs in the same region that we see a CLIP peak, we can use PyTHEA to see how, at least based on the model, it varies binding of that specific RNA binding protein by measuring the posterior of the model in the reference versus when we mutate the input. And I did this for all the polymorphisms reported in the GWAS catalog from EBI, and all of these polymorphisms are are basically annotated to be important to a specific phenotype or disease. And in some cases, I noticed there is this large deviance in the posterior of the motif. So we have some instances where a single nucleotide change can reduce the posterior of the model from 0.9 to 0.1. And these are the two most common basically reference an altered allele that is seen in the population, and the other nucleotides that are not seen in the population in this specific case, for example, they were actually not making that much of a change. So SNRpbA specifically was quite interesting because when I was quantifying the different polymorphisms according to the predicted effect that they have, I saw this quite high enrichment in splicing. Splicing polymorphisms. And SNRPB is in fact actually a splicing factor. So for this specific case, SNRPB, I wanted to see if this is like something that you would expect by chance. If you mutate like this 256-mer, would you expect to observe the same type of change or not? So I generated this null set based on mutating the same set of peaks that overlap GWAS hits and quantified how they changed the posterior of the model. And here the RBP set, which is on the y-axis, it's the real GWAS hits. You actually see they have much higher quantiles than the null, meaning there's actually an enrichment in polymorphisms, and many of these effects seem to be quite important. And here, by discussing this, I want to talk about that something like PyTHEA can have many applications in hypothesis generation. And it's not specific to just mutations in cancer or in complex disease, but it can also be applied to infectious disease. So we already know that RNA-binding proteins are quite important in the life cycle of many viruses. So here, this is a figure from the Lee et al. paper showing the polyovirus cycle and the involvement of host RNA binding proteins. And there's this complex interaction happening between the host RNA binding proteins, the viral RNA binding proteins, and proteases that are important for the viral RNA translation, viral RNA switch from translation to replication, replication in general. So I think with PyTHEA and models like PyTHEA, we can generate many hypotheses in many cases to discover more about the biology of RNA binding proteins. With that, I want to summarize that basically with PyTHEA, with introducing this fixed dilated convolution, we can capture this context-free grammar of RNA structure. And predict binding of RNA-binding proteins with good accuracy and precision. And we can use this for generating hypotheses in many cases. I would like to thank my supervisor, Hany, and also members of the Guderzi lab, and also my co-supervisor, Dr. Bo Wang, and his lab at Vector Institute. I will be happy to take any questions.