mapped_BAM_to_fastq.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Extraction of mapped reads from bam files."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```Inputs: bam files generated by STAR or other aligners```\n",
    "\n",
    "Navigate to the directory containging the BAM files.\n",
    "\n",
    "If the BAM files do now follow this naming convention: ```'*.fastq.sort.bam'```, specifier your own file identifier in the second code bock\n",
    "\n",
    "```Outputs: FASTQ files containing ONLY reads mapped to your reference genome (specified during reads alignment)```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!pwd"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "picard_jar_loc = '/home/chihchung/Installation/picard.jar'\n",
    "BAM_file_postfix = '.fastq.sort.bam'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "bamfs = [fs for fs in os.listdir('./') if fs.endswith(BAM_file_postfix) and not fs.startswith('mapped_')];bamfs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import multiprocessing\n",
    "import subprocess\n",
    "\n",
    "def work(cmd):\n",
    "    return subprocess.call(cmd, shell=True)\n",
    "count = multiprocessing.cpu_count()\n",
    "pool = multiprocessing.Pool(processes=count)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "get mapped bam files"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "if not os.path.isdir('mapped_fastq'): os.mkdir('mapped_fastq')\n",
    "cmd_get_mapped = ['samtools view -b -F 4 %s > mapped_%s'% (bam_fn, bam_fn) for bam_fn in bamfs]\n",
    "cmd_get_fastq_from_bams = [\n",
    "'java -jar %s SamToFastq \\\n",
    "     I= mapped_%s \\\n",
    "     FASTQ= mapped_fastq/%s.fastq'%(picard_jar_loc, bam_fn, bam_fn[:-len(BAM_file_postfix)]) for bam_fn in bamfs]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "pool.map(work, cmd_get_mapped)\n",
    "pool.map(work, cmd_get_fastq_from_bams)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python [conda env:pytorch]",
   "language": "python",
   "name": "conda-env-pytorch-py"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}