{'fraction_non_alphanumeric': 0.022982885085574573, 'fraction_numerical': 0.038141809290953545, 'mean_word_length': 5.032869785082174, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Epigenome-based splicing prediction using a recurrent neural network

Alternative RNA splicing provides an important means to expand metazoan transcriptome diversity. Contrary to what was accepted previously, splicing is now thought to predominantly take place during transcription. Motivated by emerging data showing the physical proximity of the spliceosome to Pol II, we surveyed the effect of epigenetic context on co-transcriptional splicing. In particular, we observed that splicing factors were not necessarily enriched at exon junctions and that most epigenetic signatures had a distinctly asymmetric profile around known splice sites. Given this, we tried to build an interpretable model that mimics the physical layout of splicing regulation where the chromatin context progressively changes as the Pol II moves along the guide DNA. We used a recurrent-neural-network architecture to predict the inclusion of a spliced exon based on adjacent epigenetic signals, and we showed that distinct spatio-temporal features of these signals were key determinants of model outcome, in addition to the actual nucleotide sequence of the guide DNA strand. After the model had been trained and tested (with >80% precision-recall curve metric), we explored the derived weights of the latent factors, finding they highlight the importance of the asymmetric time-direction of chromatin context during transcription.


Abstract 23
Alternative RNA splicing provides an important means to expand metazoan transcriptome 24 diversity. Contrary to what was accepted previously, splicing is now thought to predominantly 25 take place during transcription. Motivated by emerging data showing the physical proximity of 26 the spliceosome to Pol II, we surveyed the effect of epigenetic context on co-transcriptional 27 splicing. In particular, we observed that splicing factors were not necessarily enriched at exon 28 junctions and that most epigenetic signatures had a distinctly asymmetric profile around known 29 splice sites. Given this, we tried to build an interpretable model that mimics the physical layout 30 of splicing regulation where the chromatin context progressively changes as the Pol II moves 31 along the guide DNA. We used a recurrent-neural-network architecture to predict the inclusion 32 of a spliced exon based on adjacent epigenetic signals, and we showed that distinct spatio-33 temporal features of these signals were key determinants of model outcome, in addition to the 34 actual nucleotide sequence of the guide DNA strand. After the model had been trained and tested 35 (with >80% precision-recall curve metric), we explored the derived weights of the latent factors, 36 finding they highlight the importance of the asymmetric time-direction of chromatin context 37 during transcription. 38 39

Author Summary 40
In humans, only about 2% of the genome is comprised of so-called coding regions and can give 41 rise to protein products. However, the human transcriptome is much more diverse than the 42 number of genes found in these coding regions. Each gene can give rise to multiple transcripts 43 through a process during transcription called alternative splicing. There is a limited 44 understanding of the regulation of splicing and the underlying splicing code that determines cell-type-specific splicing. Here, we studied epigenetic features that characterize splicing regulation 46 in humans using a recurrent neural network model. Unlike feedforward neural networks, this 47 suggest that the splicing code may be directly or indirectly encoded within the epigenomic 160 context. 161 162

Clustering biosamples based on splicing patterns 163
Previous studies have shown that various epigenomic marks are correlated across similar tissues 164 and cell types [41]. It is now widely accepted that the transcriptional regulatory circuitry of a 165 particular cell type is reflected in its epigenetic landscape. To explore the potential linkage 166 between epigenetic regulation and tissue-specific splicing, we examined splicing patterns across 167 49 ENCODE biosamples. Based on a similarity of percent-splice-in (PSI) values for all coding 168 exons (n=185,405), we clustered biosamples into five categories using hierarchical clustering 169 ( Figure 3D). Splicing patterns were highly correlated among tissue types from the same cell-of-170 origin, reproducing similar clustering results based on epigenetic marks. For example, blood-171 lineage cell types formed cluster C2 whereas brain and neural cells were clustered in cluster C4. 172 Moreover, we observed that cancerous cell lines cluster together in cluster C3. 173 174 In addition to using the PSI similarity matrix to cluster cell types into categories, we can project 175 the cells onto a low-dimensional cell space using principal component analysis (PCA). We 176 measured alternative splicing patterns in terms of exonic expression level (fragment per kilobase 177 per million reads mapped, FPKM) across diverse ENCODE cell types and examined how cells 178 are placed in the context of others. Interestingly, we observed that cancer-related cell lines were 179 located proximal to each other in the PCA cell space (Supplementary Figure 4). 180 181

Modeling splicing regulation: key characteristics of an RNN architecture
To investigate the latent representation of splicing instruction encoded within the epigenomic 183 context, we aimed to construct a predictive model of splicing. We opted for an RNN architecture,184 which has proven successful in various sequential information processing and prediction tasks 185 such as natural language processing and translation [42][43][44], to explore the contribution of the 186 epigenomic context to the regulation of alternative splicing. 187 We start by describing a simple RNN, which shares many of the features we intend to model. In particular, we adopted the long short-term memory (LSTM) [45] model to describe an RNN 207 architecture. In principal, a simple RNN allows us to model a time-dependent task from 208 sequential data. However, in practice, the simple model suffers from the problem of vanishing 209 gradients, where the gradients responsible for updating weights with respect to the partial 210 derivative of error function becomes negligible in a long sequence and hampers the model from 211 learning long-term time dependencies. Therefore, we used both LSTM and gated recurrent unit 212 (GRU), which have many of the same simple intuitive properties of the simple RNN but allow 213 learning from longer sequences. The LSTM is an extension of the same idea that includes more 214 sophisticated gates, which allows the cell to retain long-term memory between cells while 215 avoiding the problem of vanishing gradients when training the network. The specific equations 216 for the LSTM model we adopted is shown in the Methods. 217 218 Modeling splicing regulation: How the RNN architecture fits the problem 219 The rationale for applying an RNN to our model is that (1) an RNN is optimized for processing 220 sequential information like genomic sequences and epigenomic profiles along genomic 221 coordinates, (2) an RNN has a time-direction resembling how RNA is transcribed by RNA  Figure 4B). The average F1 score marginally increased to 0.8573 when using 17 histone, 246 chromatin accessibility, DNA methylation, and nucleosome density profiles.
We performed the splicing prediction with or without the RBP profile and measured how much 248 predictive performance is gained from additional information. We observed a marginal 249 improvement in predictive performance when RBP binding profiles were added to the baseline 250 model (measured in improvement of F1 score from 0.84 to 0.86) (Supplementary Figure 9A, B). 251 This suggests RBP binding information may be redundant and already represented in the 252 epigenetic features. We also compared prediction results from normal cell types to those from 253 cancerous cell lines. Since previous studies on cancer-specific alternative splicing [50,51] have 254 suggested potential linkage of aberrant splicing events to the disease risk [52-55], we expected 255 to see differences in splicing regulation between normal and cancerous cell types. However, we 256 did not observe a significant difference in prediction performance between normal and cancerous 257 cell types (average F1 score for normal biosamples: 0.8465, cancerous biosamples: 0.8765). We 258 also cross-tested a model trained from one cell type to another. After we fit our model to one cell 259 type, we transferred the fitted weights and model parameters to predict splicing on other cell 260 types. When we tested between cell types from the same cell-of-origin (e.g., train on adult liver 261 model and test on HepG2 data, train on lung model and test on A549 data), we did not observe a 262 significant difference in predictive performance. However, we observed a moderate reduction in 263 splicing prediction performance when we cross-tested cells from different cell-of-origin 264 (Supplementary Figure 5B, F1 score is better metric for comparing cross-cell testing due to class 265 imbalance across cell types). Thus, the epigenomic regulatory landscape around SSs appears to 266 be generally conserved across cell types. Moreover, we compared the classification performance 267 to other models based on random forest and k-nearest neighbors and found that our model was 268 superior in terms of classification accuracy ( Figure Figure 6C). We observed that the combination of H3K36me3 and 285 H3K27ac features together contributed more than when they were used individually 286 (Supplementary Figure 6D) Figure 8D). We also noticed that H3K27me3 and K3K9me3 were negatively 306 weighted at the input gate, suggesting that these features have a negative impact on exon 307 inclusion, consistent with our previous findings. 308 309

Influence of temporal epigenetic context on splicing regulation 310
We specifically designed our splicing model to represent the physical layout of splicing 311 regulation, where a sequence of chromatin contexts is fed progressively to the model. Therefore, 312 the model takes into account the temporal direction (progression from 5' to 3' in direction). To 313 show that model has learned this asymmetric temporal relationship of epigenetic features, we first trained a baseline model (in the normal 5' to 3' direction) and then fed a series of epigenetic 315 signals in a "reverse" order (3' to 5' in direction) as input to it. We examined how the model 316 prediction behaved in this context. If the model was agnostic to the temporal direction of features, 317 both forward and reverse input features should give the same predictive power. By using a model 318

GRU model 430
We adopted the following equations for the modeling of splicing using GRU.
Pre-processing of data for the training model 436 We selected six normal and three cancer samples from the reference epigenome series. The 437 dataset contains consolidated epigenomes from the Roadmap Epigenomics Consortium [41] and 438 the ENCODE Consortium. All datasets were uniformly processed and mapped to the GRCh38 439 human reference genome. All samples contained a core set of histone modification tracks 440 (H3K4me1, H3K4me3, H3K36me3, H3K27ac, H3K27me3, and H3K9me3) as well as RNA-seq 441 data. We used additional histone modification tracks, as well as DNase I hypersensitivity, DNA 442 methylation, and nucleosome positioning tracks, to predict alternative splicing upon availability. 443 Detailed information on datasets used can be found in Supplementary Table 1 There is no single metric that can give you a measure of performance in a binary classification 461 problem. Relying on one metric can be misleading especially when there is high class imbalance. 462 Therefore, we employed various metrics to measure the performance of the predictive model. 463 ROC curve explains the tradeoff between true-positive rate (TPR) and false-positive rate (FPR). 464 PR curve visualizes the tradeoff between positive predictive value (PPV) and true-positive rate 465 (TPR). 466 In addition, we used F1-score, which is the harmonic mean of precision and recall, to measure 469 the performance of the splicing model. 470 Hyperparameter tuning of splicing model and training 472 We tested a range of dimensions and depths of RNN models and network design 473 hyperparameters to optimize the alternative splicing model. We chose optimal hyperparameters 474 by tuning one parameter at a time while fixing the rest. Hyperparameters included but were not 475 limited to the number of recurrent layers, size of neurons in each layer, pooling strategy, dropout 476 rate, choice of activation function and loss function, optimizer, and number of the epoch. We shuffled the order of the data and split the dataset into training and test sets using an 80 to 20% 478 ratio. 20% of test data was set aside for the performance evaluation. 80% of training data was 479 split again between 80 to 20% (64 and 16% of the original data) for fitting the model and 480 validating the model fit during the training phase. We fed a range of sequences from 50 to 1,000 481 bp within each SS and found the 400 bp span to be the ideal size for the model. For the neural 482 network architecture, we achieved the best result when two RNN units were stacked together, 483 which allowed the model to learn higher-level temporal representations. We used a hidden state 484 size of two by default and we recommend not using a hidden state size greater than 128 to avoid 485 overfitting problems (Supplementary Figure 8A). We applied three variants of the RNN model, We found that both LSTM and GRU were capable of learning long-term dependencies and were 489 effective in learning high-dimensional contextual relationships between epigenomic features 490 around the SSs. We split the input sequences into two parts where the first half represented a 3' 491 acceptor SS and the latter half represented a 5' donor SS. We fed these sequences into two 492 separate RNN units of size 200 and merged them into another RNN unit of size 400. The last 493 RNN layer was followed by a dropout layer to prevent overfitting of the training dataset.   3'

Chromatin Marks
Nascent RNA