File: PileupTestUtils.py

package info (click to toggle)
python-pysam 0.15.4+ds-3
  • links: PTS, VCS
  • area: main
  • in suites: bullseye, sid
  • size: 27,992 kB
  • sloc: ansic: 140,738; python: 7,881; sh: 265; makefile: 223; perl: 41
file content (160 lines) | stat: -rw-r--r-- 5,793 bytes parent folder | download | duplicates (3)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
import os
import subprocess
import pysam

from TestUtils import BAM_DATADIR, force_str

def build_pileup_with_samtoolsshell(fn):
    os.system("samtools mpileup {} 2> /dev/null | wc -l > /dev/null".format(fn))
    return 2998


def build_pileup_with_samtoolspipe(fn):
    FNULL = open(os.devnull, 'w')
    with subprocess.Popen(["samtools", "mpileup", fn],
                          stdin=subprocess.PIPE,
                          stdout=subprocess.PIPE,
                          stderr=FNULL) as proc:
        return len(proc.stdout.readlines())

    
def build_pileup_with_pysam(*args, **kwargs):
    with pysam.AlignmentFile(*args, **kwargs) as inf:
        return len(list(inf.pileup(stepper="samtools")))


def build_depth_with_samtoolsshell(fn):
    os.system(
        "samtools mpileup {} 2> /dev/null | awk '{{a += $4}} END {{print a}}' > /dev/null".format(fn))
    return 107241


def build_depth_with_samtoolspipe(fn):
    FNULL = open(os.devnull, 'w')
    with subprocess.Popen(["samtools", "mpileup", fn],
                          stdin=subprocess.PIPE,
                          stdout=subprocess.PIPE,
                          stderr=FNULL) as proc:
        data = [x.split() for x in proc.stdout.readlines()]
        return [int(x[3]) for x in data]


def build_depth_with_filter_with_pysam(*args, **kwargs):
    with pysam.AlignmentFile(*args, **kwargs) as inf:
        return [x.get_num_aligned() for x in inf.pileup(stepper="samtools")]


def build_depth_with_pysam(*args, **kwargs):
    with pysam.AlignmentFile(*args, **kwargs) as inf:
        return [x.nsegments for x in inf.pileup(stepper="samtools")]


def build_query_bases_with_samtoolsshell(fn):
    os.system("samtools mpileup {} 2> /dev/null | awk '{{a = a $5}} END {{print a}}' | wc -c > /dev/null".format(fn))
    return 116308


def build_query_bases_with_samtoolspipe(fn, *args, **kwargs):
    FNULL = open(os.devnull, 'w')
    with subprocess.Popen(["samtools", "mpileup", fn] + list(args),
                          stdin=subprocess.PIPE,
                          stdout=subprocess.PIPE,
                          stderr=FNULL) as proc:
        stdout = proc.stdout.read().decode()
        return [x.split()[4] for x in stdout.splitlines()]


def build_query_bases_with_samtoolspysam(fn, *args):
    return [x.split()[4] for x in pysam.samtools.mpileup(fn, *args).splitlines()]
    

def build_query_bases_with_pysam_pileups(*args, **kwargs):
    total_pileup = []
    with pysam.AlignmentFile(*args, **kwargs) as inf:
        total_pileup = [
            [r.alignment.query_sequence[r.query_position_or_next]
                for r in column.pileups if r.query_position_or_next is not None]
            for column in inf.pileup(stepper="samtools")]
    return total_pileup


def build_query_qualities_with_pysam_pileups(*args, **kwargs):
    total_pileup = []
    with pysam.AlignmentFile(*args, **kwargs) as inf:
        total_pileup = [
            [r.alignment.query_qualities[r.query_position_or_next]
                for r in column.pileups if r.query_position_or_next is not None]
            for column in inf.pileup(stepper="samtools")]
    return total_pileup


def build_query_bases_with_pysam(fn, *args, **kwargs):
    total_pileup = []
    with pysam.AlignmentFile(fn) as inf:
        total_pileup = [column.get_query_sequences(
            mark_ends=True, add_indels=True, mark_matches=True) for column in
                        inf.pileup(*args, **kwargs)]
    return total_pileup


def build_query_names_with_pysam(*args, **kwargs):
    total_pileup = []
    with pysam.AlignmentFile(*args, **kwargs) as inf:
        total_pileup = [column.get_query_names() for column in
                        inf.pileup(stepper="samtools")]
    return total_pileup


def build_query_qualities_with_samtoolspipe(fn):
    FNULL = open(os.devnull, 'w')
    with subprocess.Popen(["samtools", "mpileup", fn],
                          stdin=subprocess.PIPE,
                          stdout=subprocess.PIPE,
                          stderr=FNULL) as proc:
        data = [force_str(x).split()[5] for x in proc.stdout.readlines()]
    return data


def build_query_qualities_with_pysam(*args, **kwargs):
    total_pileup = []
    with pysam.AlignmentFile(*args, **kwargs) as inf:
        total_pileup = [column.get_query_qualities() for column in
                        inf.pileup(stepper="samtools")]
    return total_pileup


def build_mapping_qualities_with_samtoolspipe(fn):
    FNULL = open(os.devnull, 'w')
    with subprocess.Popen(["samtools", "mpileup", "-s", fn],
                          stdin=subprocess.PIPE,
                          stdout=subprocess.PIPE,
                          stderr=FNULL) as proc:
        data = [force_str(x).split()[6] for x in proc.stdout.readlines()]
    return data


def build_mapping_qualities_with_pysam(*args, **kwargs):
    total_pileup = []
    with pysam.AlignmentFile(*args, **kwargs) as inf:
        total_pileup = [column.get_mapping_qualities() for column in
                        inf.pileup(stepper="samtools")]
    return total_pileup


def build_query_positions_with_samtoolspipe(fn):
    FNULL = open(os.devnull, 'w')
    with subprocess.Popen(["samtools", "mpileup", "-O", fn],
                          stdin=subprocess.PIPE,
                          stdout=subprocess.PIPE,
                          stderr=FNULL) as proc:
        data = [list(map(int, force_str(x).split()[6].split(",")))
                for x in proc.stdout.readlines()]
    return data


def build_query_positions_with_pysam(*args, **kwargs):
    total_pileup = []
    with pysam.AlignmentFile(*args, **kwargs) as inf:
        total_pileup = [column.get_query_positions() for column in
                        inf.pileup(stepper="samtools")]
    return total_pileup