{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Mapping and SNP calling for all individuals"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Using a loop for mapping, with a single folder per sample"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for i in {1..20}\n",
    "    do\n",
    "        cd /ifb/data/mydatalocal/MappingAndSNP\n",
    "        echo Clone$i;\n",
    "        mkdir dirClone$i; \n",
    "        cd dirClone$i; \n",
    "        echo \"Cloning Clone$i data\";\n",
    "        ln -s /ifb/data/mydatalocal/shortReads/Clone${i}_pairedFq1.fq /ifb/data/mydatalocal/MappingAndSNP/dirClone$i/Clone${i}_1.fastq;\n",
    "        ln -s /ifb/data/mydatalocal/shortReads/Clone${i}_pairedFq2.fq /ifb/data/mydatalocal/MappingAndSNP/dirClone$i/Clone${i}_2.fastq;\n",
    "        echo -e \"\\nMapping Clone$i\\n\";\n",
    "        bwa mem -t 4 -R \"@RG\\tID:Clone${i}\\tSM:Clone${i}\" ../referenceCorrect.fasta Clone${i}_1.fastq Clone${i}_2.fastq > Clone$i.sam;\n",
    "        samtools view -@4 -b -f 0x02 Clone$i.sam | samtools sort -@4 -o Clone$i.SORTED.bam;\n",
    "        samtools index Clone$i.SORTED.bam;\n",
    "        echo -e \"\\nCalling Clone$i\";\n",
    "        gatk MarkDuplicates -I Clone$i.SORTED.bam -M duplicates.$i.metrics -O Clone$i.SORTED.MD.bam;\n",
    "        samtools index Clone$i.SORTED.MD.bam;\n",
    "        gatk --java-options \"-Xmx4g\" HaplotypeCaller --native-pair-hmm-threads 4 -I Clone$i.SORTED.MD.bam -O Clone$i.g.vcf -R ../referenceCorrect.fasta -ERC GVCF;\n",
    "    done"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Calling all samples on one raw VCF"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "cd /ifb/data/mydatalocal/MappingAndSNP\n",
    "#Loop to inflate the --variant option\n",
    "OPTION=\"\"\n",
    "for i in {1..20}\n",
    "do\n",
    "    OPTION=\"${OPTION} --variant dirClone${i}/Clone${i}.g.vcf\"\n",
    "done\n",
    "echo $OPTION\n",
    "gatk CombineGVCFs -R referenceCorrect.fasta $OPTION -O rawSNP.vcf"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Have a Look to it"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "head -n 1000 rawSNP.vcf | tail"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Compute the Genotypes *:D I know I am a bad person...*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gatk --java-options \"-Xmx4g\" GenotypeGVCFs -R referenceCorrect.fasta -V rawSNP.vcf -O output.vcf"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Compute the SNP density along the chromosomes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "echo -e \"Reference\\t1000000\\n\" > /ifb/data/mydatalocal/MappingAndSNP/genome.txt\n",
    "bedtools genomecov -bga -split -i /ifb/data/mydatalocal/MappingAndSNP/output.vcf -g /ifb/data/mydatalocal/MappingAndSNP/genome.txt > /ifb/data/mydatalocal/MappingAndSNP/density.csv"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "head /ifb/data/mydatalocal/MappingAndSNP/density.csv"
   ]
  }
 ],
 "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": 4
}