{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "5e58d68b-9b88-4446-ada8-04271bf4d511",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    },
    "tags": []
   },
   "source": [
    "## Réaliser la partie bioinformatique d'une analyse RNAseq\n",
    "\n",
    "### Récupération des données RNAseq :\n",
    "Les fichiers fastq ont été préparés (téléchargés depuis l'ENA et réduits) et sont déposés sur les ressources du cluster de l'IFB, dans l'espace projet de la formation. \n",
    "\n",
    "Comme nous venons de copier l'ensemble du dossier dans votre espace projet, les données ont été copiées et sont dans le dossier `atelier_rnaseq/01-Bioinfo/data/`.\n",
    "\n",
    "Vérifiez que vous êtes bien dans votre répertoire copié et sinon, déplacez-vous y."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "13f4ae88-4807-4e5a-88ac-f1424bb6a13e",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    },
    "tags": []
   },
   "outputs": [],
   "source": [
    "# verifiez que vous etes dans votre copie:\n",
    "pwd"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0d70d7e2-0c92-4414-a623-a92d93c3297a",
   "metadata": {},
   "outputs": [],
   "source": [
    "# sinon, allez-y en remplacant <votre_espace_projet> \n",
    "# (il ne doit plus y avoir de \"<\" ni \">\")\n",
    "cd /shared/projects/<votre_espace_projet>/atelier_rnaseq/01-Bioinfo"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4dbca273-2d2e-4bca-a1e0-4ed386f7605d",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    },
    "tags": []
   },
   "source": [
    "Listez les données RNAseq, il y a 3 fichiers WT et 3 KO au format fastq:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d48d4bf6-a61a-45f4-9d08-fcce36ede83d",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    },
    "tags": []
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "0ac031dd-78df-4645-9864-f599b6e9c989",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    },
    "tags": []
   },
   "source": [
    "Créez un dossier, results, pour placer les résultats"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4ba8bdbc-27de-40eb-9c67-ec275977cacd",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    },
    "tags": []
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "738ee09a-3f8d-40e0-9f3d-1b4c7923bd12",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    },
    "tags": []
   },
   "source": [
    "Les différentes étapes d'une analyse RNAseq sont :\n",
    "- analyse qualité des séquences\n",
    "- nettoyage (si besoin)\n",
    "- mapping des reads sur le génome\n",
    "- comptage des reads sur chaque gènes\n",
    "\n",
    "La première étape d'analyse de la qualité des séquences se lance avec le logiciel FastQC.\n",
    "Nous appliquons le nettoyage des reads (même si pas besoin, juste pour montrer) avec le logiciel Fastp.\n",
    "Pour le mapping, nous utilisons l'outil STAR, version 2.7.9a.\n",
    "Pour changer de format, sélectionner des alignements, ou indexer le fichier bam de mapping pour le visualisé, nous utilisons samtools (version 1.15.1) qui permet de travailler avec les fichiers résultants du mapping.\n",
    "Pour le comptage, nous utilisons featureCounts pour le comptage (featureCounts appartient à un ensemble d'outils nommé subread, nous utilisons la version 1.6.1) \n",
    "\n",
    "Sur le cluster de calcul de l'IFB, il faut activer les outils dont nous avons besoin.\n",
    "Ainsi, nous avons 5 outils à activer avec la commande \"module\" :\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9c90e6f3-29a0-4f02-83cd-a99fa110d162",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    },
    "tags": []
   },
   "outputs": [],
   "source": [
    "## Charger les outils dont tu as besoin pour ton analyse\n",
    "module load fastqc/0.11.9 #Controle qualite des sequences\n",
    "module load fastp/0.23.1 #Nettoyage des sequences (si besoin)\n",
    "module load star/2.7.9a #Mapping des reads sur le genome\n",
    "module load samtools/1.15.1 #Selection des alignements\n",
    "module load subread/1.6.1 #Comptage avec featureCounts"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "26c185a1-d46a-4ab4-af10-fd318e00b907",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    },
    "tags": []
   },
   "source": [
    "La documentation pour fastQC est là: https://www.bioinformatics.babraham.ac.uk/projects/fastqc/.\n",
    "\n",
    "On va lancer fastQC sur un échantillon (KO1) en utilisant le répertoire results:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8480e708-5f55-4c58-af07-e6ac76ba0cbe",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    },
    "tags": []
   },
   "outputs": [],
   "source": [
    "fastqc -h #pour voir toutes les options"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5406ad5a-5a5e-47a0-bab1-7144d2bd7430",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    },
    "tags": []
   },
   "source": [
    "On vérifie les fichiers résultats et en particulier le fichier .html qui va vous donner un recap des résultats.\n",
    "Si besoin (ici pas vraiment besoin), on passe au nettoyage:"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3d6867cc-5b51-44dd-8a9d-547d1466eb4f",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    },
    "tags": []
   },
   "source": [
    "La documentation pour fastp est ici: https://github.com/OpenGene/fastp\n",
    "ou:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "dc6ea8a7-5c49-43aa-8e75-38995eff3864",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    },
    "tags": []
   },
   "outputs": [],
   "source": [
    "fastp -h #pour voir toutes les options"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d3be379b-64e5-43ca-b82b-81f9647c0507",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    },
    "tags": []
   },
   "source": [
    "Pour lancer sur l'échantillon KO1:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "be318a1e-8814-4c86-903f-86061a0c1a0e",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    },
    "tags": []
   },
   "outputs": [],
   "source": [
    "fastp -i data/KO1_R1.fastq.gz -I data/KO1_R2.fastq.gz \\\n",
    "      -o results/KO1_R1.clean.fastq.gz -O results/KO1_R2.clean.fastq.gz \\\n",
    "      --detect_adapter_for_pe \n",
    "#Et d'autres paramètres si besoin en fonction du résultat du fastqc"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "deaad819-56d3-4e8d-8fab-5e1d6929cf26",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    },
    "tags": []
   },
   "source": [
    "On revérifie les séquences nettoyées avec fastqc"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cb179e33-dfb3-402f-a452-babff9870ab9",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    },
    "tags": []
   },
   "outputs": [],
   "source": [
    "fastqc -h #pour voir toutes les options\n",
    "fastqc -o results results/KO1_R1.clean.fastq.gz\n",
    "fastqc -o results results/KO1_R2.clean.fastq.gz"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3b057b34-dd35-425e-a28e-1bad7d5d77a2",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    },
    "tags": []
   },
   "source": [
    "On peut maintenant passer au mapping des reads sur le génome.\n",
    "La documentation de l'outil STAR est ici https://physiology.med.cornell.edu/faculty/skrabanek/lab/angsd/lecture_notes/STARmanual.pdf\n",
    "\n",
    "L'usage de base de STAR est constituée de 2 étapes :\n",
    "- indexation du génome\n",
    "- mapping des reads\n",
    "\n",
    "L'indexation est un procédé de compression et de ré-organisation des séquences qui permet aux logiciels de gagner beaucoup de temps. Chaque outil a le sien (blast aussi).\n",
    "\n",
    "Le génome d'Arabidopsis thaliana a déjà été indexé et est accessible dans l'espace des banques de séquence : /shared/bank/arabidopsis_thaliana/TAIR10.1/star-2.7.9a/\n",
    "\n",
    "Si votre génome n'est pas disponible, l'option \"--runMode genomeGenerate\" de STAR permet de réaliser l'indexation du génome pour STAR."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c0cec75f-0074-4bda-8cac-1a65da4c2732",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    },
    "tags": []
   },
   "source": [
    "Les options de bases pour le lancement de STAR sont décrites ici :\n",
    "https://physiology.med.cornell.edu/faculty/skrabanek/lab/angsd/lecture_notes/STARmanual.pdf#subsection.3.1\n",
    "\n",
    "- NumberOfThreads : --runThreadN entier\n",
    "- accès au génome indexé : --genomeDir /path/to/genomeDir\n",
    "- fichier de reads : --readFilesIn /path/to/read1 [/path/to/read2 ]\n",
    "- pour des fichiers compressés (.gz) : --readFilesCommand zcat \n",
    "\n",
    "Nous devons donc les mentionner et adapter leur valeur dans la ligne de lancement.\n",
    "\n",
    "STAR engendre beaucoup de fichiers en sortie, et l'option STAR --outFileNamePrefix /path/to/output/dir/prefix permet de directmenent les ranger dans un sous-dossier. Avec \"prefix\" il est possible de spécifier quelques lettres qui correpondent à l'échantillon. \n",
    "Par ex. pour l'échantillon WT1, on peut écrire results/WT1_ \n",
    "\n",
    "STAR a été entraîné par défaut pour le génome humain. Pour _Arabidopsis thaliana_, nous pouvons (devons) modifier les valeurs par défaut de certaines options : \n",
    "- maximum intron size --alignIntronMax 1000 (la distribution de la longueur moyenne des introns des gènes d'Arabidopsis thaliana est bimodale, à 80 et à 120 bp), mieux que ~600 Kb par défaut\n",
    "- maximum gap between two mates --alignMatesGapMax 10000 (defaut ~600 Kb)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2bd97dd6-409b-40e0-b2b9-74bf51f8c708",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    },
    "tags": []
   },
   "source": [
    "Voilà, nous avons tout réuni pour écrire (puis lancer) la commande de mapping. \n",
    "Continuons avec l'échantillon \"KO1\" :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7f961991-313e-45f8-8e54-8e53a00c6f26",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    },
    "tags": []
   },
   "outputs": [],
   "source": [
    "STAR --genomeDir /shared/bank/arabidopsis_thaliana/TAIR10.1/star-2.7.9a/ --runThreadN 4 \\\n",
    "     --readFilesCommand zcat --outFileNamePrefix results/KO1_ \\\n",
    "     --readFilesIn results/KO1_R1.clean.fastq.gz results/KO1_R2.clean.fastq.gz \\\n",
    "     --outSAMtype BAM SortedByCoordinate --alignIntronMax 1000 --alignMatesGapMax 10000 \\\n",
    "     --sjdbGTFfile /shared/bank/arabidopsis_thaliana/TAIR10.1/gtf/GCF_000001735.4_TAIR10.1_genomic.gtf\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c73774ad-e058-469f-9438-2f9c9aa1936b",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    },
    "tags": []
   },
   "source": [
    "Une fois la commande exécutée, il y a plusieurs fichiers dans le dossier \"results\". Le fichier de sortie d'alignement est celui qui se termine par .bam."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2bab39ae-11ea-40d7-a4c5-e1c7b2c5b568",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    },
    "tags": []
   },
   "source": [
    "Ce fichier bam n'est pas indexé pour être correctement lu par IGV. Pour cela on utilise samtools index (documentation de samtools : https://arq5x.github.io/tutorials/samtools/samtools.html#samtools-index):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3a108041-a2ea-40d9-9d06-a729fffc1f79",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    },
    "tags": []
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "83030a93-a947-4b1a-b206-8c4727400dd1",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    },
    "tags": []
   },
   "source": [
    "Et enfin, on peut passer au comptage du nombre de reads alignés sur chaque gene avec featureCounts. Chercher dans la documentation (https://subread.sourceforge.net/featureCounts.html) les options pour :\n",
    "- compter sur les exons\n",
    "- afficher le nom des gènes\n",
    "- compter en mode orienté\n",
    "- utiliser 4 threads"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "97d669ab-7f10-416d-856d-180c03826fdf",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    },
    "tags": []
   },
   "outputs": [],
   "source": [
    "featureCounts ..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c827feb7-5917-40d0-804a-3b67db2a107b",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    },
    "tags": []
   },
   "source": [
    "Si on récapitule, voici toutes les commandes à lancer pour un échantillon:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3d622d51-32dc-42aa-9b54-8c12adc2b1a4",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    },
    "tags": []
   },
   "outputs": [],
   "source": [
    "#Controle qualite\n",
    "fastqc -o results data/KO1_R1.fastq.gz\n",
    "fastqc -o results data/KO1_R2.fastq.gz\n",
    "#Nettoyage\n",
    "fastp -i data/KO1_R1.fastq.gz -I data/KO1_R2.fastq.gz \\\n",
    "      -o results/KO1_R1.clean.fastq.gz -O results/KO1_R2.clean.fastq.gz \\\n",
    "      --detect_adapter_for_pe #Et d'autres paramètres si besoin en fonction du résultat du fastqc\n",
    "#Controle qualite sur fichiers nettoyes\n",
    "fastqc -o results results/KO1_R1.clean.fastq.gz\n",
    "fastqc -o results results/KO1_R2.clean.fastq.gz\n",
    "#Mapping\n",
    "STAR --genomeDir /shared/bank/arabidopsis_thaliana/TAIR10.1/star-2.7.9a/ --runThreadN 4 \\\n",
    "     --readFilesCommand zcat --outFileNamePrefix results/KO1_ \\\n",
    "     --readFilesIn results/KO1_R1.clean.fastq.gz results/KO1_R2.clean.fastq.gz \\\n",
    "     --outSAMtype BAM SortedByCoordinate --alignIntronMax 1000 --alignMatesGapMax 10000 \\\n",
    "     --sjdbGTFfile /shared/bank/arabidopsis_thaliana/TAIR10.1/gtf/GCF_000001735.4_TAIR10.1_genomic.gtf\n",
    "#indexation\n",
    "samtools index results/KO1_Aligned.sortedByCoord.out.bam\n",
    "#Comptage des reads\n",
    "featureCounts -T 4 -t exon -g gene_id -s 1 \\\n",
    "              -a /shared/bank/arabidopsis_thaliana/TAIR10.1/gtf/GCF_000001735.4_TAIR10.1_genomic.gtf \\\n",
    "              -o results/KO1_feature.out results/KO1_Aligned.sortedByCoord.out.bam 2> results/KO1_counting.logs\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5e2cbcff-f233-411d-bd66-12b99ed7f758",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    },
    "tags": []
   },
   "source": [
    "Allons-y pour 2 KO et 2 WT :D, ça donne envie?? (Ne lancez pas ça, on vous montrera plus tard comment faire beaucoup plus simple :))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f8338f72-7160-4e5e-a3c7-9697700dcb24",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    },
    "tags": []
   },
   "outputs": [],
   "source": [
    "#KO1\n",
    "fastqc -o results data/KO1_R1.fastq.gz\n",
    "fastqc -o results data/KO1_R2.fastq.gz\n",
    "fastp -i data/KO1_R1.fastq.gz -I data/KO1_R2.fastq.gz \\\n",
    "      -o results/KO1_R1.clean.fastq.gz -O results/KO1_R2.clean.fastq.gz \\\n",
    "      --detect_adapter_for_pe \n",
    "fastqc -o results results/KO1_R1.clean.fastq.gz\n",
    "fastqc -o results results/KO1_R2.clean.fastq.gz\n",
    "STAR --genomeDir /shared/bank/arabidopsis_thaliana/TAIR10.1/star-2.7.9a/ --runThreadN 4 \\\n",
    "     --readFilesCommand zcat --outFileNamePrefix results/KO1_ \\\n",
    "     --readFilesIn results/KO1_R1.clean.fastq.gz results/KO1_R2.clean.fastq.gz \\\n",
    "     --outSAMtype BAM SortedByCoordinate --alignIntronMax 1000 --alignMatesGapMax 10000 \\\n",
    "     --sjdbGTFfile /shared/bank/arabidopsis_thaliana/TAIR10.1/gtf/GCF_000001735.4_TAIR10.1_genomic.gtf\n",
    "samtools index results/KO1_Aligned.sortedByCoord.out.bam\n",
    "featureCounts -T 4 -t exon -g gene_id -s 1 \\\n",
    "              -a /shared/bank/arabidopsis_thaliana/TAIR10.1/gtf/GCF_000001735.4_TAIR10.1_genomic.gtf \\\n",
    "              -o results/KO1_feature.out results/KO1_Aligned.sortedByCoord.out.bam 2> results/KO1_counting.logs\n",
    "\n",
    "#KO2\n",
    "fastqc -o results data/KO2_R1.fastq.gz\n",
    "fastqc -o results data/KO2_R2.fastq.gz\n",
    "fastp -i data/KO2_R1.fastq.gz -I data/KO2_R2.fastq.gz \\\n",
    "      -o results/KO2_R1.clean.fastq.gz -O results/KO2_R2.clean.fastq.gz \\\n",
    "      --detect_adapter_for_pe \n",
    "fastqc -o results results/KO2_R1.clean.fastq.gz\n",
    "fastqc -o results results/KO2_R2.clean.fastq.gz\n",
    "STAR --genomeDir /shared/bank/arabidopsis_thaliana/TAIR10.1/star-2.7.9a/ --runThreadN 4 \\\n",
    "     --readFilesCommand zcat --outFileNamePrefix results/KO2_ \\\n",
    "     --readFilesIn results/KO2_R1.clean.fastq.gz results/KO2_R2.clean.fastq.gz \\\n",
    "     --outSAMtype BAM SortedByCoordinate --alignIntronMax 1000 --alignMatesGapMax 10000 \\\n",
    "     --sjdbGTFfile /shared/bank/arabidopsis_thaliana/TAIR10.1/gtf/GCF_000001735.4_TAIR10.1_genomic.gtf\n",
    "samtools index results/KO2_Aligned.sortedByCoord.out.bam\n",
    "featureCounts -T 4 -t exon -g gene_id -s 1 \\\n",
    "              -a /shared/bank/arabidopsis_thaliana/TAIR10.1/gtf/GCF_000001735.4_TAIR10.1_genomic.gtf \\\n",
    "              -o results/KO2_feature.out results/KO2_Aligned.sortedByCoord.out.bam 2> results/KO2_counting.logs\n",
    "\n",
    "\n",
    "#WT1\n",
    "fastqc -o results data/WT1_R1.fastq.gz\n",
    "fastqc -o results data/WT1_R2.fastq.gz\n",
    "fastp -i data/WT1_R1.fastq.gz -I data/WT1_R2.fastq.gz \\\n",
    "      -o results/WT1_R1.clean.fastq.gz -O results/WT1_R2.clean.fastq.gz \\\n",
    "      --detect_adapter_for_pe \n",
    "fastqc -o results results/WT1_R1.clean.fastq.gz\n",
    "fastqc -o results results/WT1_R2.clean.fastq.gz\n",
    "STAR --genomeDir /shared/bank/arabidopsis_thaliana/TAIR10.1/star-2.7.9a/ --runThreadN 4 \\\n",
    "     --readFilesCommand zcat --outFileNamePrefix results/WT1_ \\\n",
    "     --readFilesIn results/WT1_R1.clean.fastq.gz results/WT1_R2.clean.fastq.gz \\\n",
    "     --outSAMtype BAM SortedByCoordinate --alignIntronMax 1000 --alignMatesGapMax 10000 \\\n",
    "     --sjdbGTFfile /shared/bank/arabidopsis_thaliana/TAIR10.1/gtf/GCF_000001735.4_TAIR10.1_genomic.gtf\n",
    "samtools index results/WT1_Aligned.sortedByCoord.out.bam\n",
    "featureCounts -T 4 -t exon -g gene_id -s 1 \\\n",
    "              -a /shared/bank/arabidopsis_thaliana/TAIR10.1/gtf/GCF_000001735.4_TAIR10.1_genomic.gtf \\\n",
    "              -o results/WT1_feature.out results/WT1_Aligned.sortedByCoord.out.bam 2> results/WT1_counting.logs\n",
    "\n",
    "#KO2\n",
    "fastqc -o results data/WT2_R1.fastq.gz\n",
    "fastqc -o results data/WT2_R2.fastq.gz\n",
    "fastp -i data/WT2_R1.fastq.gz -I data/WT2_R2.fastq.gz \\\n",
    "      -o results/WT2_R1.clean.fastq.gz -O results/WT2_R2.clean.fastq.gz \\\n",
    "      --detect_adapter_for_pe \n",
    "fastqc -o results results/WT2_R1.clean.fastq.gz\n",
    "fastqc -o results results/WT2_R2.clean.fastq.gz\n",
    "STAR --genomeDir /shared/bank/arabidopsis_thaliana/TAIR10.1/star-2.7.9a/ --runThreadN 4 \\\n",
    "     --readFilesCommand zcat --outFileNamePrefix results/WT2_ \\\n",
    "     --readFilesIn results/WT2_R1.clean.fastq.gz results/WT2_R2.clean.fastq.gz \\\n",
    "     --outSAMtype BAM SortedByCoordinate --alignIntronMax 1000 --alignMatesGapMax 10000 \\\n",
    "     --sjdbGTFfile /shared/bank/arabidopsis_thaliana/TAIR10.1/gtf/GCF_000001735.4_TAIR10.1_genomic.gtf\n",
    "samtools index results/WT2_Aligned.sortedByCoord.out.bam\n",
    "featureCounts -T 4 -t exon -g gene_id -s 1 \\\n",
    "              -a /shared/bank/arabidopsis_thaliana/TAIR10.1/gtf/GCF_000001735.4_TAIR10.1_genomic.gtf \\\n",
    "              -o results/WT2_feature.out results/WT2_Aligned.sortedByCoord.out.bam 2> results/WT2_counting.logs\n",
    "\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Bash",
   "language": "bash",
   "name": "bash"
  },
  "language_info": {
   "codemirror_mode": "shell",
   "file_extension": ".sh",
   "mimetype": "text/x-sh",
   "name": "bash"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}