sv_length_stats.py

import argparse

import matplotlib.pyplot as plt
import numpy as np

from general_functions import load_sv_conflict
from sv_type_stat import get_sv_type

def parse_args():
    parser = argparse.ArgumentParser(formatter_class=argparse.ArgumentDefaultsHelpFormatter)
    parser.add_argument("conflicts", type=str)
    parser.add_argument("--sv_type", type=str, default=None)
    parser.add_argument("--file", type=str, default=None, help="If given, puts the conflicts with chosen sv_type and avg length lower than upper_bound into new file")
    parser.add_argument("--fig", type=str, default=None, help="If given, saves a histogramm of the distribution of avg lengths to the file.")
    parser.add_argument("--upper_bound", "-u", type=int, help="Only conflicts with an avg length below this value are considered", default=10000)
    parser.add_argument("--lower_bound", "-l", type=int, help="Only conflicts with an avg length above this value are considered", default=50)
    parser.add_argument("--bins", "-b", type=int, default=20, help="Bins for the histogramm.")

    args = parser.parse_args()
    return args


if __name__ == "__main__":
    args = parse_args()
    conflict_file = open(args.conflicts)

    avg_lengths = []
    if args.file:
        out_file = open(args.file, "w")

    for conflict in conflict_file:
        svs = load_sv_conflict(conflict)

        lengths = []
        for sv in svs:
            if args.sv_type and not get_sv_type(sv) == args.sv_type:
                break
            sv_length = abs(int(sv.split('_')[2]))
            lengths.append(sv_length)
        else:
            avg_len = sum(lengths) / len(lengths)
            avg_lengths.append(avg_len)
            if args.file and args.lower_bound < avg_len < args.upper_bound:
                out_file.write(conflict)

    if args.fig:
        logbins = np.logspace(np.log10(args.lower_bound), np.log10(args.upper_bound), args.bins)
        plt.hist(avg_lengths, bins=logbins)
        plt.xscale('log')
        plt.xlabel("Avg length of calls per conflict")
        plt.ylabel("Frequency")
        plt.savefig(args.fig)