{"id":737,"date":"2019-06-05T15:24:54","date_gmt":"2019-06-05T13:24:54","guid":{"rendered":"https:\/\/blogs.fu-berlin.de\/gruenstaeudl\/?p=737"},"modified":"2019-06-05T15:29:28","modified_gmt":"2019-06-05T13:29:28","slug":"extracting-mapped-reads-from-sam-file","status":"publish","type":"post","link":"https:\/\/blogs.fu-berlin.de\/gruenstaeudl\/2019\/06\/05\/extracting-mapped-reads-from-sam-file\/","title":{"rendered":"Extracting mapped R1 and R2 reads from SAM file"},"content":{"rendered":"<p><strong>Extracting those that mapped<\/strong><\/p>\n<p>Recently, I found myself in need of extracting only those reads from a sequence alignment map (SAM) file that actually mapped to the reference genome, while maintaining the separation into the paired-end read design. By using a combination of <a href=\"http:\/\/samtools.sourceforge.net\/\" target=\"_blank\" rel=\"noopener\">samtools<\/a>, <a href=\"https:\/\/bedtools.readthedocs.io\/en\/latest\/\" target=\"_blank\" rel=\"noopener\">bedtools<\/a> and <em>awk<\/em>, this can be done very efficiently:<\/p>\n<div style=\"background-color: #ffebdb\">\n<pre>INFL=MySamples.sam\r\nSTEM=${INF%.sam*}\r\nTMP1=${STEM}.extracted.sam\r\nTMP2=${STEM}.sorted.sam\r\nOTFL=${STEM}.fastq\r\n\r\n# Extracting only paired mapped reads\r\nsamtools view -b -F12 $INFL &gt; $TMP1\r\n#samtools view -b -F4 $INFL &gt; $TMP1\r\n\r\n# Sorting SAM file by read header\r\n# (necessary for the command `bedtools bamtofastq`)\r\nsamtools sort -n $TMP1 &gt; $TMP2\r\n\r\n# Extracting fastq sequences\r\nbedtools bamtofastq -i $TMP2 -fq $OTFL\r\n\r\n# Separating successfully mapped paired reads \r\n# into an R1 file and an R2 file\r\ncat $OTFL | awk -v n=4 'BEGIN{f=2} { \r\n  if((NR-1)%n==0){f=1-f}; print &gt; \"out_R\" f \".fastq\"\r\n}'\r\nmv out_R-1.fastq ${STEM}_R1.fastq\r\nmv out_R2.fastq ${STEM}_R2.fastq\r\n\r\n<\/pre>\n<\/div>\n<p>&nbsp;<\/p>\n","protected":false},"excerpt":{"rendered":"<p>Extracting those that mapped Recently, I found myself in need of extracting only those reads from a sequence alignment map (SAM) file that actually mapped to the reference genome, while maintaining the separation into the paired-end read design. By using a combination of samtools, bedtools and awk, this can be done very efficiently: INFL=MySamples.sam STEM=${INF%.sam*} [&hellip;]<\/p>\n","protected":false},"author":2306,"featured_media":0,"comment_status":"open","ping_status":"open","sticky":false,"template":"","format":"standard","meta":{"footnotes":""},"categories":[57598],"tags":[],"class_list":["post-737","post","type-post","status-publish","format-standard","hentry","category-bioinformatics"],"_links":{"self":[{"href":"https:\/\/blogs.fu-berlin.de\/gruenstaeudl\/wp-json\/wp\/v2\/posts\/737","targetHints":{"allow":["GET"]}}],"collection":[{"href":"https:\/\/blogs.fu-berlin.de\/gruenstaeudl\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/blogs.fu-berlin.de\/gruenstaeudl\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/blogs.fu-berlin.de\/gruenstaeudl\/wp-json\/wp\/v2\/users\/2306"}],"replies":[{"embeddable":true,"href":"https:\/\/blogs.fu-berlin.de\/gruenstaeudl\/wp-json\/wp\/v2\/comments?post=737"}],"version-history":[{"count":6,"href":"https:\/\/blogs.fu-berlin.de\/gruenstaeudl\/wp-json\/wp\/v2\/posts\/737\/revisions"}],"predecessor-version":[{"id":743,"href":"https:\/\/blogs.fu-berlin.de\/gruenstaeudl\/wp-json\/wp\/v2\/posts\/737\/revisions\/743"}],"wp:attachment":[{"href":"https:\/\/blogs.fu-berlin.de\/gruenstaeudl\/wp-json\/wp\/v2\/media?parent=737"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/blogs.fu-berlin.de\/gruenstaeudl\/wp-json\/wp\/v2\/categories?post=737"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/blogs.fu-berlin.de\/gruenstaeudl\/wp-json\/wp\/v2\/tags?post=737"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}