Identification d’un épissage alternatif et d’une polyadénylation dans les données de séquençage de l’ARN

Ce protocole fournit une compréhension complète des isoformes géniques générées par l’épissage alternatif et la polyadénylation en fournissant un flux de travail étape par étape pour identifier les sites d’épissage différentiel, les exons exprimés différemment et les sites poly(A). Le principal avantage de ce protocole est qu’il évalue à la fois les méthodes basées sur les exons et les méthodes basées sur les événements pour étudier l’épissage alternatif. Il applique également une méthode à base d’exons pour étudier la polyadénylation alternative.

Les fichiers R Markdown qui incluent les codes et les notes pour l’analyse AS et AP ont été fournis. Il serait conseillé de suivre les étapes du fichier R Markdown et d’atteindre attentivement la note pour chaque étape. Pour identifier l’épissage différentiel à l’aide de diffSplice from limma, suivez le fichier R notebook.

Préparez les fichiers d’entrée comme décrit dans le manuscrit texte. Assurez-vous que les étapes un à trois du manuscrit ont été suivies séquentiellement pour préparer les fichiers d’entrée avant de continuer. Commencez par charger les bibliothèques nécessaires.

Pour effectuer un filtrage non spécifique, extrayez d’abord la matrice des comptes de lecture obtenus précédemment et créez une liste d’entités à l’aide de la fonction DGEList à partir du package edgeR, où les lignes représentent les gènes et les colonnes représentent les échantillons. Ensuite, transformez les données de l’échelle brute en nombre par million à l’aide de la fonction CPM du package edgeR et conservez les exons dont le nombre est supérieur à un seuil définissable. Cet ensemble de données contient six échantillons.

Par conséquent, la MPC est fixée à plus d’un échantillon et au moins trois échantillons sur six. Normalisez les comptages entre les échantillons à l’aide de la fonction calcNormFactors du package edgeR à l’aide de la moyenne tronquée des valeurs M. Cette fonction calcule les facteurs d’échelle pour ajuster la taille des bibliothèques.

Utilisez la table d’exemples générée précédemment pour créer la matrice de plan afin de définir les conditions expérimentales pour chaque échantillon. Exécutez la fonction voom du package limma pour traiter les données de séquençage de l’ARN afin d’estimer la variance. Cette fonction générera des poids de précision pour corriger le bruit de comptage de Poisson et transformer le nombre de niveaux d’exons en nombre de deux comptes par million ou logCPM.

Exécutez la fonction lmfit pour ajuster les modèles linéaires aux données d’expression de chaque exon. Exécutez ensuite la fonction eBayes pour calculer des statistiques empiriques pour le modèle ajusté afin de détecter l’expression différentielle des exons. Définir une matrice de contraste pour les comparaisons expérimentales d’intérêt.

Utilisez les contrastes. Fonction d’ajustement pour obtenir des coefficients et des erreurs-types pour chaque paire de comparaisons. Exécutez diffSplice sur le modèle ajusté pour tester les différences dans l’utilisation des gènes par exon entre le type sauvage et le knockout.

Explorez les résultats les mieux classés à l’aide de la fonction topSplice où un test égal à t donne un classement des exons AS et un test égal à simes donne un classement des gènes. Exécutez la fonction plotSplice pour tracer les résultats. En mettant le gène d’intérêt dans l’argument de l’identification du gène, les points rouges montrent les exons exprimés différemment.

Générez un tracé de volcan à l’aide du package bioconducteur EnhancedVolcano pour afficher les exons exprimés différentiellement. Pour utiliser rMATS, assurez-vous que la dernière version de rMATS version 4.1.1 est installée à l’aide de conda ou GitHub dans le répertoire de travail. Accédez au dossier contenant les fichiers bam obtenus après mappage.

Préparez les fichiers texte requis par rMATS pour les deux conditions de copie du nom des fichiers bam et de leur chemin séparés par une virgule. Exécutez rmas. PY à l’aide des deux fichiers texte d’entrée générés décrivant le chemin des fichiers BAM et l’annotation.

GTF obtenu précédemment. Cela génère un dossier de sortie contenant rmats_out fichiers texte décrivant les statistiques, y compris les valeurs P et les niveaux d’inclusion pour chaque événement d’épissage séparément. Utilisez le maser du boîtier bioconducteur pour explorer les résultats rMATS.

Chargez les fichiers texte de jonction et de comptage d’exons avec l’extension JCEC dans l’objet maser et incluez au moins cinq lectures moyennes par événement d’épissage pour filtrer le résultat en fonction de la couverture. Pour visualiser les résultats rMATS, exécutez d’abord la fonction topEvents à partir du package maser, en sélectionnant les événements d’épissage significatifs à un taux de fausse découverte de 10 % et une variation minimale de 10 % du pourcentage d’épissage dans ou PSI. Vérifiez les événements génétiques pour les gènes individuels d’intérêt et tracez les valeurs PSI pour chaque événement d’épissage de ce gène.

Générez un tracé de volcan en spécifiant le type d’événement. Utilisez les résultats des événements d’épissage obtenus avec rMATS sous forme de fichiers texte pour générer des tracés sashimi à l’aide du package rmats2sashimiplot. Le diagramme de sashimi montre un événement d’exon sauté dans le gène Wnk1.

Chaque rangée représente un échantillon de séquençage d’ARN, trois répétitions de type sauvage et un knockout Mbnl1. La hauteur montre la couverture de lecture en RPKM et les arcs de connexion représentent les lectures de jonction à travers les exons. La partie inférieure montre des isoformes alternatives annotées du modèle de gènes.

Un changement de pli substantiel et de fortes preuves statistiques de différences réelles peuvent être observés dans les gènes situés dans les quadrants supérieurs gauche ou droit des parcelles volcaniques obtenues en utilisant diffSplice et DEXSeq. Un exon de cassette variait entre différentes conditions pour le gène Wnk1. Le graphique d’utilisation différentielle des exons a montré des signes d’épissage différentiel sur cinq sites d’exons près de Wnk1.6.45, les exons surlignés en rose étant susceptibles d’être épissés dans les échantillons knockout Mbnl1 par rapport au type sauvage.

Le tracé volcanique des gènes qui sont alternativement épissés a aidé à distinguer les gènes qui ont été exclus du type sauvage et ceux qui ont été inclus dans le type sauvage. Les types d’épissage SE, A5SS, A3SS, MXE et RI ont été visualisés à l’aide de diagrammes de sashimi des gènes les plus significatifs de ces événements. L’activité différentielle de l’APA dans trois régions principales non traduites de gènes a été observée à l’aide de diagrammes volcaniques.

Les résultats d’utilisation du site PA significativement différentiels acquis à partir de différents pipelines ont été visualisés à l’aide d’un diagramme d’événements. Un déplacement distal à proximal significatif de l’utilisation du site PA dans les doubles knockouts peut être observé dans les deux gènes FOSL1 et Papola. La couverture moyenne dans les régions flanquantes ancrées à des sites de clivage PA connus à l’échelle du génome a été déterminée à l’aide d’un diagramme diagnostique.

Assurez-vous que les paramètres tels que les informations transspécifiques et l’autorisation du chevauchement multiple sont correctement utilisés lors de la génération de mesures de comptage. L’ajustement du modèle linéaire et la génération de paires de contraste sont importants pour une comparaison correcte. Pour rMATS, assurez-vous que tous les paramètres sont correctement définis en fonction de vos données avant d’exécuter la commande.

Les gènes obtenus à partir de l’activité d’épissage différentiel pourraient être utilisés pour effectuer une analyse d’enrichissement de l’ensemble de gènes. Un autre outil appelé MISO pourrait être utilisé pour une analyse plus approfondie basée sur les événements.