<commit_before>"""
===============================================================================
Use pysam model to calculate the AS from bam file
===============================================================================
Author : <NAME>
Date   : 2014-03-25 14:29:08
"""
import sys
import re
import optparse
import os
import string
import pysam 
import matplotlib.pyplot as plt
import numpy as np
import AlterAlign as ATA

def IsSNP(refbase, alleles):

    isSnp = True
    for ale in alleles: 
        if len(ale) > 1 or len(refbase) > 1: isSnp = False

    return isSnp

def main(opt): 

    vcfInfile = opt.vcfInfile
    bamInfile = opt.bamInfile
    faInfile  = opt.refInfile
    sampleID  = opt.sample
    refId     = opt.refChrId
    outPrefix = opt.outPrefix
    mapq      = opt.mapq
    newID     = opt.exc
    fa = ATA.LoadFaSeq(faInfile)

    print >> sys.stderr, '# [INFO] Fa Loading finish ***'
    if bamInfile[-4:] == '.bam':
        samInHandle = pysam.Samfile(bamInfile, 'rb')
    else :
        samInHandle = pysam.Samfile(bamInfile, 'r')

    #samOutHandle = pysam.Samfile(outPrefix + '.bam', 'wb', template=samInHandle)
    vcfOutHandle = open(outPrefix + '.vcf', 'w')

    print >> sys.stderr, '# [INFO] Now Scaning the VCF and doing alternate align ... ...'
    if vcfInfile[-3:] == '.gz':
        I = os.popen('gzip -dc %s' % vcfInfile)
    else :
        I = open(vcfInfile)

    frist = True
    while 1:
        lines = I.readlines(100000)
        if not lines : break;

        for line in lines :
        
            line = line.strip('\n')
            col  = line.split()
            if re.search(r'^##', line): 
                if frist and re.search(r'^##FORMAT=', line): 
                    frist = False
                    vcfOutHandle.write('##FORMAT=<ID=AA,Number=4,Type=Integer,Description="Information of Alternate Align. Format: Ref_perfect,Alt_Perfect,Both_Perfect,Both_Imperfect">\n')
                vcfOutHandle.write('%s\n' % line)
                continue
            elif re.search(r'^#CHROM', line):<commit_msg>add tabix to speed up the process in reading vcf file with specific chromoID<commit_after>"""
===============================================================================
Use pysam model to calculate the AS from bam file
===============================================================================
Author : <NAME>
Date   : 2014-03-25 14:29:08
"""
import sys
import re
import optparse
import os
import string
import pysam 
import matplotlib.pyplot as plt
import numpy as np
import AlterAlign as ATA

def IsSNP(refbase, alleles):

    isSnp = True
    for ale in alleles: 
        if len(ale) > 1 or len(refbase) > 1: isSnp = False

    return isSnp

def main(opt): 

    vcfInfile = opt.vcfInfile
    bamInfile = opt.bamInfile
    faInfile  = opt.refInfile
    sampleID  = opt.sample
    refId     = opt.refChrId
    outPrefix = opt.outPrefix
    mapq      = opt.mapq
    newID     = opt.exc
    fa = ATA.LoadFaSeq(faInfile, refId)

    print >> sys.stderr, '# [INFO] Fa Loading finish ***'
    if bamInfile[-4:] == '.bam':
        samInHandle = pysam.Samfile(bamInfile, 'rb')
    else :
        samInHandle = pysam.Samfile(bamInfile, 'r')

    #samOutHandle = pysam.Samfile(outPrefix + '.bam', 'wb', template=samInHandle)
    vcfOutHandle = open(outPrefix + '.vcf', 'w')

    print >> sys.stderr, '# [INFO] Now Scaning the VCF and doing alternate align ... ...'
    if vcfInfile[-3:] == '.gz':
        if refId == "ALL":
            I = os.popen('gzip -dc %s' % vcfInfile)
        else:
            I = os.popen('/home/siyang/Bin/software_pip/tabix-0.2.6/tabix -h %s %s' % (vcfInfile, refId))
    else :
        I = open(vcfInfile)

    frist = True
    while 1:
        lines = I.readlines(100000)
        if not lines : break;

        for line in lines :
        
            line = line.strip('\n')
            col  = line.split()
            if re.search(r'^##', line): 
                if frist and re.search(r'^##FORMAT=', line): 
                    frist = False
                    vcfOutHandle.write('##FORMAT=<ID=AA,Number=4,Type=Integer,Description="Information of Alternate Align. Format: Ref_perfect,Alt_Perfect,Both_Perfect,Both_Imperfect">\n')
                vcfOutHandle.write('%s\n' % line)
                continue
            elif re.search(r'^#CHROM', line):
{'fraction_non_alphanumeric': 0.18860555938433265, 'fraction_numerical': 0.012634964392373075, 'mean_word_length': 2.511290322580645, 'pattern_counts': {'":': 1, '<': 7, '<?xml version=': 0, '>': 19, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, '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': '2482515', 'n_tokens_mistral': 1499, 'n_tokens_neox': 1381, 'n_words': 356}