<commit_before>"""



http://www.acgt.me/blog/2014/12/16/understanding-mapq-scores-in-sam-files-does-37-42

http://biofinysics.blogspot.fr/2014/05/how-does-bowtie2-assign-mapq-scores.html


https://gitlab.univ-nantes.fr/a-slide/ContaVect/blob/9a411abfa720064c205c5f6c811afdfea206ed12/pyDNA/pySamTools/Bam.py

pysamtools
"""

import pandas as pd
import pylab
import pysam

# pysam uses htlib behing the scene and is very fast
# pysam works great for BAM file but with SAM, it needs to read the file after
# each compete iteration, which is not very useful


class BAM(pysam.AlignmentFile):
    """

    mode rb for bam files

    """
    def __init__(self, filename, mode="rb", *args):
        # works for py27 but not py35 probably a missing __init__  or __new__ in
        # AlignmentFile class. See e.g., stackoverflow/questions/
        # 26653401/typeerror-object-takes-no-parameters-but-only-in-python-3
        #super(BAM, self).__init__(filename, mode, *args)
        pysam.AlignmentFile.__init__(filename, mode, *args)

    def get_read_names(self):
        self.reset()
        names = [this.qname for this in self]
        self.reset()
        return names

    def iter_unmapped_reads(self):
        self.reset()
        unmapped = (this.qname for this in self if this.is_unmapped)
        self.reset()
        return unmapped

    def iter_mapped_reads(self):
        self.reset()
        mapped = (this.qname for this in self if this.is_unmapped is False)
        self.reset()
        return mapped

    def __len__(self):
        self.reset()
        N = len([x for x in self])
        self.reset()
        return N

    def get_stats(self):
        d = {}
        d['total_reads'] = len(list(self.iter_unmapped_reads()))
        d['mapped_reads'] = len(list(self.iter_mapped_reads()))
        d['unmapped_reads'] = len(list(self.iter_unmapped_reads()))
        d['contamination [%]'] = float(d['mapped_reads']) /float(d['unmapped_reads']) 
        d['contamination [%]'] *= 100
        return d

class SAMFlags(object):
<commit_msg>add get_flags utilities<commit_after>"""



http://www.acgt.me/blog/2014/12/16/understanding-mapq-scores-in-sam-files-does-37-42

http://biofinysics.blogspot.fr/2014/05/how-does-bowtie2-assign-mapq-scores.html


https://gitlab.univ-nantes.fr/a-slide/ContaVect/blob/9a411abfa720064c205c5f6c811afdfea206ed12/pyDNA/pySamTools/Bam.py

pysamtools



Interesting commands::

    samtools flagstat contaminant.bam
    samtools stats contaminant.bam
"""

import pandas as pd
import pylab
import pysam

# pysam uses htlib behing the scene and is very fast
# pysam works great for BAM file but with SAM, it needs to read the file after
# each compete iteration, which is not very useful


class BAM(pysam.AlignmentFile):
    """

    mode rb for bam files

    """
    def __init__(self, filename, mode="rb", *args):
        # works for py27 but not py35 probably a missing __init__  or __new__ in
        # AlignmentFile class. See e.g., stackoverflow/questions/
        # 26653401/typeerror-object-takes-no-parameters-but-only-in-python-3
        #super(BAM, self).__init__(filename, mode, *args)
        pysam.AlignmentFile.__init__(filename, mode, *args)

    def get_read_names(self):
        self.reset()
        names = [this.qname for this in self]
        self.reset()
        return names

    def iter_unmapped_reads(self):
        self.reset()
        unmapped = (this.qname for this in self if this.is_unmapped)
        self.reset()
        return unmapped

    def iter_mapped_reads(self):
        self.reset()
        mapped = (this.qname for this in self if this.is_unmapped is False)
        self.reset()
        return mapped

    def __len__(self):
        self.reset()
        N = len([x for x in self])
        self.reset()
        return N

    def get_stats(self):
        #See samtools stats

        """1526795 + 0 in total (QC-passed reads + QC-failed reads)
        13 + 0 secondary
        0 + 0 supplementary
        0 + 0 duplicates
        3010 + 0 mapped (0.20% : N/A)
        1526782 + 0 paired in sequencing
        763391 + 0 read1
        763391 + 0 read2
        2700 + 0 properly paired (0.18% : N/A)
        2976 + 0 with itself and mate mapped
        21 + 0 singletons (0.00% : N/A)
        0 + 0 with mate mapped to a different chr
        0 + 0 with mate mapped to a different chr (mapQ>=5)
        """

        d = {}
        d['total_reads'] = len(list(self.iter_unmapped_reads()))
        d['mapped_reads'] = len(list(self.iter_mapped_reads()))
        d['unmapped_reads'] = len(list(self.iter_unmapped_reads()))
        d['contamination [%]'] = float(d['mapped_reads']) /float(d['unmapped_reads']) 
        d['contamination [%]'] *= 100
        return d

    def get_flags(self):
        self.reset()
        flags = [s.flag for s in self]
        self.reset()
        return flags

    def get_flags_as_df(self):
        """

            b = BAM(filename)
            df = b.get_flags_as_df()
            df.sum()

        """
        flags = self.get_flags()
        data = [(this, [flag&this for flag in flags]) 
            for this in [1,2,4,8,16,32,64,128,256,512,1024,2048]]
        import pandas as pd
        df = pd.DataFrame(dict(data))
        df = df > 0
        return df



class SAMFlags(object):

{'fraction_non_alphanumeric': 0.11534831036435718, 'fraction_numerical': 0.04058901264866906, 'mean_word_length': 2.429126213592233, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 5, 'https://': 2, 'lorem ipsum': 0, 'www.': 2, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '2087424', 'n_tokens_mistral': 1966, 'n_tokens_neox': 1762, 'n_words': 440}