NNTM_accuracy.sagews

len('a8babb7a161b06722aa93a585031d1703ecad697')

40

print 1+1
reset()
print 1+2

2 3

from rna_poly import RNAPolytope
from glob import glob
from os.path import basename,splitext
from datetime import datetime
import pickle
import __main__
from basic_sec_struct import sec_struct
__main__.RNAPolytope = RNAPolytope
load('rna_poly.py')
load('polytope_plots.sage')
load('main.sage')

get_bname = lambda fname:splitext(basename(fname))[0]
def timestamp(message=''):
    return "[{0}] {1}".format(datetime.now(),message)
def load_obj(name):
    with open(name, 'rb') as f:
        return pickle.load(f)
def my_hot(x):
    return colormaps['hot_r'](x/2+0.25)
def my_Greys(x):
    return colormaps['Greys'](x/4+0.5)
pkl_names = glob.glob('/home/user/data/bio/tRNA/st_tuples/*.pkl')

names=[ get_bname(n)[:-3] for n in pkl_names]
sobj_names = ['/home/user/data/bio/tRNA/sobjs/{0}.sobj'.format(n) for n in names]
#names

from itertools import combinations,chain
from numpy import average as avg, median as med, std
from collections import defaultdict

def get_sobj_name(name):
    return '/home/user/data/bio/tRNA/sobjs/{0}.sobj'.format(name)

def get_suboptimal_ranking(name,P,b):
    D = load_obj("/home/user/data/bio/tRNA/st_tuples/"+name+"_ss.pkl")
    K=D.keys()
    region_density = {}
    slices = sliced_cones_b(P, b)
    relevant_signatures = [tuple(sl.original_vertex) for sl in slices]
    max_relevant_subopt = 0
    total_relevant_structs=0
    for k in K:
        #k is a tuple: (a,b,c,d,default_ss)
        #D[k] is a list of (ss,(x,y,z,w,t1,t2,t3,t4,acc))'s
        for S,(x,y,z,w,t1,t2,t3,t4,acc) in D[k]:
            break
        x,y,z,w = (QQ(e) for e in (x,y,z,w))
        if (x,y,z,w) in relevant_signatures:
            n=len(D[k])
            total_relevant_structs +=n
            region_density[(x,y,z,w)] = n
            if n > max_relevant_subopt:
                max_relevant_subopt = n

    for k in region_density:
        region_density[k] /= float(max_relevant_subopt)

    total_structures = sum([len(D[k]) for k in K])
    max_structures = max([len(D[k]) for k in K])
    return region_density, (total_relevant_structs,total_structures), (max_relevant_subopt,max_structures)

def get_average_accuracy_ranking(name):
    D = load_obj("/home/user/data/bio/tRNA/st_tuples/"+name+"_ss.pkl")
    K=D.keys()
    avg_accuracy = {}
    max_accuracy = {}
    b_avg = {}
    b_max = {}
    for k in K:
        #k is a tuple: (a,b,c,d,default_ss)
        #D[k] is a list of (ss,(x,y,z,w,t1,t2,t3,t4,acc))'s
        avg_acc = 0
        n = len(D[k])
        max_acc = 0
        #print D[k][0]
        for S,(x,y,z,w,t1,t2,t3,t4,acc) in D[k]:
            avg_acc += acc
            if acc>max_acc:
                max_acc = acc
        avg_acc /= n
        x,y,z,w = (QQ(e) for e in (x,y,z,w))
        avg_accuracy[(x,y,z,w)] = avg_acc
        b_avg[avg_acc] = QQ(k[1])
        max_accuracy[(x,y,z,w)] = max_acc
        b_max[max_acc] = QQ(k[1])
    return avg_accuracy,b_avg

remove_dangles = lambda s: s.replace('<','.').replace('<','.')

def compute_subopt_similarities(name,omit_singletons=True):
    D = load_obj("/home/user/data/bio/tRNA/st_tuples/"+name+"_ss.pkl")
    K=D.keys()
    min_accs = {}
    signs = []
    n_subopts = []
    structs={}
    for k in K:
        subopts = list(set([remove_dangles(s[0]) for s in D[k]]))
        n_subopts.append(len(subopts))
        _,(x,y,z,w,_,_,_,_,_) = D[k][0]
        if n_subopts[-1]==1:
            if omit_singletons:
                continue
            min_acc = 1.0
        else:
            min_acc = min([min([f_measure(s,t) for s,t in combinations(subopts,2)]),min([f_measure(t,s) for s,t in combinations(subopts,2)])])
        min_accs[(x,y,z,w)]=min_acc
        signs.append((x,y,z,w))
        structs[(x,y,z,w)]=subopts
    print timestamp("{0}: max_subopts={1}, avg_subopts={2}; std_subopts={3}".format(name,max(n_subopts),avg(n_subopts),std(n_subopts)))
    return min_accs,structs,signs
def get_subopt_structures(name,default_ss,sign=False):
    D = load_obj("/home/user/data/bio/tRNA/st_tuples/"+name+"_ss.pkl")
    K=D.keys()
    for k in K:
        ss =  k[4]
        if ss==default_ss:
            if sign:
                S,(x,y,z,w,t1,t2,t3,t4,acc) = D[k][0]
                return [s[0] for s in D[k]],(x,y,z)
            else:
                return [x[0] for x in D[k]]
    #print "Returning no structures?"
    if sign:
        return [], (0,0,0)
    return []
def get_accuracies_from_struct(name,default_ss):
    D = load_obj("/home/user/data/bio/tRNA/st_tuples/"+name+"_ss.pkl")
    K=D.keys()
    for k in K:
        ss =  k[4]
        if ss==default_ss:
            return [x[1][-1] for x in D[k]]

def get_accuracies(name,default_ss):
    D = load_obj("/home/user/data/bio/tRNA/st_tuples/"+name+"_ss.pkl")
    K=D.keys()
    for k in K:
        ss =  k[4]
        if ss==default_ss:
            n = len(D[k])
            return [acc for S,(x,y,z,w,t1,t2,t3,t4,acc) in D[k]]
    return None

classic_pt = (34/10,0,4/10)
common_pt = (489/40, 51/320, -231/80)

def get_classic_slice(P,name):
    classic_slices=[]
    for sl in P.d1_slices():
        if (34/10,0,4/10) in sl:
            classic_slices.append(sl)

    return classic_slices

def classic_slices_and_accuracy(names):
    classic_slices = {}
    accuracies = {}
    for i in range(len(names)):
        name = names[i]
        sobj_name = get_sobj_name(name)

        P =RNAPolytope.construct_from_file(sobj_name)
        classic_slices[name] = get_classic_slice(P,name)
        accs = []
        for sl in classic_slices[name]:
            ss = sl.original_structure
            accs.extend(get_accuracies(name,ss))
        accuracies[name] = accs

    return classic_slices, accuracies

def get_slice_containing_pt(P,pt):
    slices = []
    for sl in P.d1_slices():
        if pt in sl:
            slices.append(sl)
    return slices

def pt_slices_and_accuracy(names,pt):
    slices={}
    accuracies={}
    for i in range(len(names)):
        name = names[i]
        sobj_name = get_sobj_name(name)

        P =RNAPolytope.construct_from_file(sobj_name)
        slices[name] = get_slice_containing_pt(P,pt)
        accs = []
        for sl in slices[name]:
            ss = sl.original_structure
            accs.extend(get_accuracies(name,ss))
        accuracies[name] = accs

    return slices, accuracies



def get_max_average_slice(P,name):
    avg_accuracy, b_avg = get_average_accuracy_ranking(name)
    max_avg_acc = 0
    max_sigs = []
    for sig in avg_accuracy:
        if avg_accuracy[sig]>max_avg_acc:
            max_avg_acc = avg_accuracy[sig]
            max_sigs = [sig]
        elif avg_accuracy[sig]==max_avg_acc:
            max_sigs.append(sig)
    out_slices = []
    for sl in P.d1_slices():
        x,y,z,w = sl.original_vertex
        if (x,y,z,w) in max_sigs:
            out_slices.append(sl)
    return out_slices

def get_neighboring_slices(P,sl):
    neighbors = []
    bounded_nbrs = []
    for s in P.d1_slices():
        if s==sl:
            print "Avoiding self"
            continue
        inter = sl.intersection(s)
        if not inter.is_empty():
            neighbors.append(s)
            if len(inter.rays())==0:
                bounded_nbrs.append(s)
    return neighbors, bounded_nbrs

def f_measure(s,t):
    S = sec_struct(s)
    T = sec_struct(t)
    acc = S.relative_accuracy(T)
    return acc

def get_similarity_matrices(P,name,sl):
    nbrs,bdd_nbrs = get_neighboring_slices(P,sl)
    subopts = {}

    for nbr in nbrs:
        subopts[nbr] = get_subopt_structures(name,nbr.original_structure)

    subopts[sl] = get_subopt_structures(name,sl.original_structure)

    matrices={}
    bdd_matrices={}
    for nbr in nbrs:
        matrix = []
        for s in subopts[nbr]:
            row = []
            for t in subopts[sl]:
                row.append(f_measure(s,t))
            matrix.append(row)
        matrices[nbr] = matrix
        if nbr in bdd_nbrs:
            bdd_matrices[nbr]=matrix
    return matrices,bdd_matrices

def get_neighbors_accuracies(P,name,sl):
    nbrs,bdd_nbrs = get_neighboring_slices(P,sl)
    accuracies={}
    bdd_accuracies={}
    for nbr in nbrs:
        accs = get_accuracies(name, nbr.original_structure)
        matrices[nbr] = accs
        if nbr in bdd_nbrs:
            bdd_accuracies[nbr]=accs
    return accuracies,bdd_accuracies

def neighbors_accuracies(names,pt):
    accuracies = {}
    bdd_accuracies = {}
    for i in range(len(names)):
        print timestamp("{0}".format(i+1)),
        name = names[i]
        sobj_name = get_sobj_name(name)

        accuracies[name]={}
        bdd_accuracies[name]={}

        P =RNAPolytope.construct_from_file(sobj_name)
        slice_list = get_slice_containing_pt(P,pt)
        for sl in slice_list:
            if len(sl.rays())==0:
                accuracies[name][sl],_ = get_neighbors_accuracies(P,name,sl)
            else:
                accuracies[name][sl],bdd_accuracies[name][sl] = get_neighbor_accuracies(P,name,sl)
    return accuracies,bdd_accuracies

def similarity_matrices(names,pt):
    matrices = {}
    bdd_matrices = {}
    for i in range(len(names)):
        print timestamp("{0}".format(i)),
        name = names[i]
        sobj_name = get_sobj_name(name)

        matrices[name]={}
        bdd_matrices[name]={}

        P =RNAPolytope.construct_from_file(sobj_name)
        slice_list = get_slice_containing_pt(P,pt)
        for sl in slice_list:
            if len(sl.rays())==0:
                matrices[name][sl],_ = get_similarity_matrices(P,name,sl)
            else:
                matrices[name][sl],bdd_matrices[name][sl] = get_similarity_matrices(P,name,sl)
    return matrices,bdd_matrices

def get_max_average_slices(names=names):
    slices = {}
    sl2na = {}
    for i in range(len(names)):
        name = names[i]
        sobj_name = get_sobj_name(name)

        P =RNAPolytope.construct_from_file(sobj_name)
        P_slices = get_max_average_slice(P,name)
        slices[name]=P_slices
        for sl in P_slices:
            sl2na[sl] = name
        if len(P_slices)>1:
            print "*"*40
            print "Found more than one max avg slice in {0}".format(name)
    return slices, sl2na

def max_avg_intersection_graph(slices):
    slice_to_label = {sl:i for i,sl in enumerate(slices)}
    label_to_slice = {i:sl for i,sl in enumerate(slices)}
    n = len(slice_to_label.keys())
    G=Graph(n)
    edges = []
    for sl1,sl2 in combinations(slices,2):
        inter = sl1.intersection(sl2)
        # Bug! Change this to inter.is_empty()
        if len(inter.rays())>0 or inter.volume()>0:
            edges.append([slice_to_label[sl1],slice_to_label[sl2]])
    G.add_edges(edges)
    return G,label_to_slice

def greedy_clique_partition(G):
    H=G.copy()
    cliques = []
    while H.num_verts()>0:
        c = H.clique_maximum()
        H.delete_vertices(c)
        cliques.append(c)
    return cliques

def plot_max_avg_acc_slice(name,P,colormap):
    avg_accuracy, b_avg = get_average_accuracy_ranking(name)
    max_avg = max(avg_accuracy.values())
    return accuracy_cone_plots_b(P,avg_accuracy,b=b_avg[max_avg],abounds=(-200,200),cbounds=(-200,200),plot_title=name,colormap=colormap)#.show(figsize=12)

def plot_suboptimal_density_slice(name,P,b=0,colormap=my_hot):
    region_density, totals, maxes = get_suboptimal_ranking(name,P,b)
    title = '{0}; total structs: {1}; max per region: {3}'.format(name, totals[0],totals[1], maxes[0],maxes[1])
    return accuracy_cone_plots_b(P,region_density,b=b,abounds=(-200,200),cbounds=(-200,200),plot_title=title,colormap=colormap)#.show(figsize=12)

def plot_max_avg_slices(names):
    for i in range(len(names)):
        name = names[i]
        sobj_name = get_sobj_name(name)

        P = P=RNAPolytope.construct_from_file(sobj_name)
        plot_max_avg_acc_slice(name,P)

def slice_center(sl,delta=1/2):
    vertex_sum = vector(sl.base_ring(), [0]*sl.ambient_dim())
    for v in sl.vertex_generator():
        vertex_sum += v.vector()
    vertex_sum=vertex_sum / (sl.n_vertices())
    for v in sl.ray_generator():
        vertex_sum += delta*v.vector()
    vertex_sum.set_immutable()
    return vertex_sum

classic = (34/10,0,4/10)
def get_classic_slices(P):
    slices = []
    for i,sl in enumerate(P.d1_slices()):
        if classic in sl:
            slices.append(sl)
    return slices

plot_slice = lambda sl: sum(sl.plot().all[1:])

def get_slice_width_height_depth_intervals(sl):
    from sage.numerical.mip import MIPSolverException
    I = []
    for i in range(3):
        lp,V=sl.to_linear_program(return_variable=True)
        var = V[i]
        lp.set_objective(var)
        try:
            o1 = lp.solve()
        except MIPSolverException as e:
            if str(e).find("unbounded")>-1:
                o1 = infinity
            else:
                raise e
        lp.set_objective(-var)
        try:
            o2 = -lp.solve()
        except MIPSolverException as e:
            if str(e).find("unbounded")>-1:
                o2 = -infinity
            else:
                print e
        I.append((o2,o1))
    return I

def is_bounded(sl):
    return len(sl.rays())==0
def is_unbounded(sl):
    return len(sl.rays())>0
def is_wedge(sl):
    return len(sl.rays())==2
def is_strip(sl):
    return len(sl.rays())==1

def get_slice_width_height_depth(sl):
    from sage.numerical.mip import MIPSolverException
    I = []
    for i in range(3):
        lp,V=sl.to_linear_program(return_variable=True)
        var = V[i]
        lp.set_objective(var)
        try:
            o1 = lp.solve()
        except MIPSolverException as e:
            if str(e).find("unbounded")>-1:
                o1 = infinity
            else:
                raise e
        lp.set_objective(-var)
        try:
            o2 = -lp.solve()
        except MIPSolverException as e:
            if str(e).find("unbounded")>-1:
                o2 = -infinity
            else:
                print e
        I.append(o1-o2)
    return I

def get_slice_width_interval(sl):
    return get_slice_width_height_depth(sl)[0]
def get_slice_height_interval(sl):
    return get_slice_width_height_depth(sl)[1]
def get_slice_depth_interval(sl):
    return get_slice_width_height_depth(sl)[2]

def get_slice_width(sl):
    a,b=get_slice_width_interval(sl)
    return b-a
def get_slice_height(sl):
    a,b=get_slice_height_interval(sl)
    return b-a
def get_slice_depth(sl):
    a,b=get_slice_depth_interval(sl)
    return b-a


def get_var_intervals(sl,pt=None):
    if pt is None:
        pt = slice_center(sl)
    elif pt=='classic':
        pt=classic
    from sage.numerical.mip import MIPSolverException
    I = []
    for i in range(3):
        lp,V=sl.to_linear_program(return_variable=True)
        var = V[i]
        lp.set_objective(var)
        for j in range(3):
            if j!=i:
                lp.add_constraint(V[j]==pt[j])
        try:
            o1 = lp.solve()
        except MIPSolverException as e:
            if str(e).find("unbounded")>-1:
                o1 = infinity
            else:
                print e
        lp.set_objective(-var)
        try:
            # Check if this should have a negative sign!
            o2 = lp.solve()
        except MIPSolverException as e:
            if str(e).find("unbounded")>-1:
                o2 = -infinity
            else:
                print e
        w = min([abs(o1-pt[i]),abs(o2-pt[i])])
        I.append(w)
    return I

def get_bounded_partitions(sname,pct=False):
    P=RNAPolytope.construct_from_file(sname)
    bounded = defaultdict(lambda:0)
    dirs = 'abc'
    for sl in P.d1_slices():
        if sl.dimension()!=3:
            continue
        rays = sl.rays()
        unb = ''
        for i in range(3):
            unb= unb + (dirs[i] if any([v[i] for v in rays]) else '')

        bounded[unb] += 1
        if unb=='a':
            print timestamp("Unbounded only in a: {0}; {1}; {2}".format(sname, sl.original_vertex, slice_center(sl)))
        elif unb=='c':
            print timestamp("Unbounded only in c: {0}; {1}; {2}".format(sname, sl.original_vertex, slice_center(sl)))
        elif unb=='ab':
            print timestamp("Unbounded only in ab: {0}; {1}; {2}".format(sname, sl.original_vertex, slice_center(sl)))
        elif unb == 'bc':
            print timestamp(sname)

    n = len(P.d1_slices())

    total_bounded = bounded['']
    total_unbounded = n-total_bounded

    if pct:
        for k in bounded:
            bounded[k] *= 100.0/n

    bounded['bounded'] = (100.0*total_bounded)/n
    bounded['unbounded'] = (100.0*total_unbounded)/n

    return bounded

P=RNAPolytope.construct_from_file(sobj_names[0])
for sl in P.d1_slices():
    break

is_bounded(sl)

True

sorted(set([1,5,7,3,5,7,3,4,7,8,3,5,7]))

[1, 3, 4, 5, 7, 8]

def timestamp_wrapper(function):
    def wrapped_function(*args, **kwargs):
        print(timestamp("[{0}] Start".format(function.__name__[:10])))
        rv = function(*args, **kwargs)
        print(timestamp("[{0}] Stop".format(function.__name__[:10])))
        return rv
    return wrapped_function

def emin(X):
    try:
        m = min(X)
    except ValueError:
        m = Infinity
    return m

def emax(X):
    try:
        M = max(X)
    except ValueError:
        M = -Infinity
    return M

def intersects_positive_orthant(sl):
    positive_orthant = Polyhedron(ieqs=[[0, 1, 0, 0], [0, 0, 1, 0], [0, 0, 0, 1]],base_ring=QQ)
    intersection = positive_orthant.intersection(sl)
    return intersection.dimension()==3

# We only compute the intersection with the vertical slab in the positive orthant where
# 0 \leq b \leq 1
def positive_orthant_intersection(sl):
    positive_orthant = Polyhedron(ieqs=[[0, 1, 0, 0], [0, 0, 1, 0], [0, 0, 0, 1],[1,0,-1,0]],base_ring=QQ)
    intersection = positive_orthant.intersection(sl)
    return intersection

def collect_positive_orthant_signatures(P):
    signature_dict = defaultdict(list)
    volume = 0
    for sl in P.d1_slices():
        intersection = positive_orthant_intersection(sl)
        if intersection.dimension()==3:
            x,y,z,w = sl.original_vertex
            signature_dict[(x,z)].append(sl)
            if is_unbounded(intersection) and (x!=0 or z!=0):
                print "Found non (0,0) unbuonded region intersecting the positive orthant!!"
                print "Signature is: {0}".format(sl.original_vertex)
                print "Rays: {0}".format(sl.rays_list())
            if is_bounded(intersection):
                volume += intersection.volume()

    return volume,signature_dict

def count_positive_orthant_slices(P):
    counter = 0
    for sl in P.d1_slices():
        if intersects_positive_orthant(sl):
            counter+= 1
    return counter

def coordinate_extremes(sl):
    vxs = sl.vertices_list()
    rays = [r for r in sl.rays_list()]

    max_a = emax([v[0] for v in vxs])
    min_a = emin([v[0] for v in vxs])
    max_b = emax([v[1] for v in vxs])
    min_b = emin([v[1] for v in vxs])
    max_c = emax([v[2] for v in vxs])
    min_c = emin([v[2] for v in vxs])

    # The following loop only affects the unbounded regions.
    for v in rays:
        a,b,c = v
        if a>0:
            max_a =  Infinity
        if a<0:
            min_a = -Infinity
        if b>0:
            max_b =  Infinity
        if b<0:
            min_b = -Infinity
        if c>0:
            max_c =  Infinity
        if c<0:
            min_c = -Infinity

    return ((min_a,max_a),(min_b,max_b),(min_c,max_c))

max_slice_norm = lambda sl: max([v.vector().norm(1) for v in sl.vertices()])
max_polytope_norm = lambda P: max([max_slice_norm(sl) for sl in P.d1_slices() if is_bounded(sl)])
max_polytope_norm_non1 = lambda P: max([max_slice_norm(sl) for sl in P.d1_slices() if is_bounded(sl) and nx1filter(sl)])
x1filter = lambda sl: sl.original_vertex.vector()[0]==1
nx1filter = lambda sl: sl.original_vertex.vector()[0]>1
true_filter = lambda x: True

@timestamp_wrapper
def get_bounded_region_intervals(P,slice_filter=true_filter):
    As=[]
    Bs=[]
    Cs=[]

    for sl in P.d1_slices():
        if slice_filter(sl):
            if is_bounded(sl):
                (a,A),(b,B),(c,C) = coordinate_extremes(sl)
                As.append((a,A))
                Bs.append((b,B))
                Cs.append((c,C))
    return As,Bs,Cs

@timestamp_wrapper
def get_bounded_region_widths(P,slice_filter=true_filter):
    As=[]
    Bs=[]
    Cs=[]

    for sl in P.d1_slices():
        if slice_filter(sl):
            if is_bounded(sl):
                (a,A),(b,B),(c,C) = coordinate_extremes(sl)
                As.append(A-a)
                Bs.append(B-b)
                Cs.append(C-c)
    return As,Bs,Cs

@timestamp_wrapper
def get_slice_change_values(P, slice_filter=true_filter):
    As = []
    Bs = []
    Cs = []

    for sl in P.d1_slices():
        if slice_filter(sl):
            (a,A),(b,B),(c,C) = coordinate_extremes(sl)
            As.extend([a,A])
            Bs.extend([b,B])
            Cs.extend([c,C])

    As = sorted(set(As))
    Bs = sorted(set(Bs))
    Cs = sorted(set(Cs))

    return As,Bs,Cs

def projected_slice(sl,b):
    b = QQ(b)
    eqn = (-b, 0, 1, 0)

    vecprojector = lambda vec: [vec[0], vec[2]]

    subspace = Polyhedron(eqns = [eqn] )

    slice2 = sl.intersection(subspace)
    if slice2.dimension() < 2:
        return None

    # The resulting coneslice is a two-dimensional polyhedron, but it still lives in Q3, so let's project down to Q2
    newverts = lambda polyslice: [vecprojector(vert) for vert in polyslice.vertices()] # Projects the defining vertices of a slice into Q2
    newrays = lambda polyslice: [vecprojector(ray) for ray in polyslice.rays()] # Projects the defining rays of a slice into Q2
    polyprojector = lambda polyslice: Polyhedron(vertices = newverts(polyslice), rays = newrays(polyslice)) # Projects a polytope from Q3 to Q2 using the helper functions we just defined

    # We can then apply polyprojector to our slice to move it down to Q2
    projslice = polyprojector(slice2)
    return projslice

from itertools import combinations
distance = lambda u,v: sum([abs(x-y) for x,y in zip(u,v)])
diameter = lambda sl: max([distance(u,v) for u,v in combinations(sl.vertices_list(),2)])


def poly_has_positive_slope(polygon):
    hull=cyclic_sort_vertices_2d(polygon.vertices())
    for i,v in enumerate(hull):
        u = hull[i-1]
        if slope(u,v)>0:
            return True
    return False
def slice_has_positive_slope(sl):
    b_vals = [p[1] for p in sl.vertices_list()]
    if is_unbounded(sl):
        try:
            diam = diameter(sl)
        except ValueError as e:
            diam = 1
            print sl.original_vertex
            print sl.original_structure
        for r in sl.ray_generator():
            v = r.vector()
            if v.norm(1)>diam:
                b_vals.append(v[1])
            else:
                w = max([1,diam])*v
                b_vals.append(w[1])
    b_coords = sorted(set(b_vals))
    n_bs = len(b_coords)
    b_cuts = [(b_coords[i]+b_coords[i-1])/2 for i in range(1,n_bs)]
    for b in b_cuts:
        poly = projected_slice(sl,b)
        if poly is not None:
            if not is_bounded(poly):
                continue
            if poly_has_positive_slope(poly):
                return True
    return False

def find_slice(P,criteria):

    for sl in P.d1_slices():
        if criteria(sl):
            return sl
    return None


def plot_slice_cuts(sl):
    b_vals = [p[1] for p in sl.vertices_list()]
    b_coords = sorted(set(b_vals))
    n_bs = len(b_coords)
    b_cuts = [(b_coords[i]+b_coords[i-1])/2 for i in range(1,n_bs)]
    for b in b_cuts:
        poly = projected_slice(sl,b)
        if poly is not None:
            poly.show()
def timestamp_wrapper(function):
    def wrapped_function(*args, **kwargs):
        print(timestamp("[{0}] Start".format(function.__name__[:10])))
        rv = function(*args, **kwargs)
        print(timestamp("[{0}] Stop".format(function.__name__[:10])))
        return rv
    return wrapped_function

@timestamp_wrapper
def collect_positive_slope_regions(P):
    n_regions = len(P.d1_slices())
    n_bounded = 0
    n_bounded_positive = 0
    n_unbounded = 0
    n_unbounded_positive = 0
    # bounded_positive_slope_slices = []
    # unbounded_positive_slope_slices = []
    for sl in P.d1_slices():
        if is_unbounded(sl):
            n_unbounded += 1
            if slice_has_positive_slope(sl):
                # unbounded_positive_slope_slices.append(sl)
                n_unbounded_positive+=1
        else:
            n_bounded += 1
            if slice_has_positive_slope(sl):
                n_bounded_positive+=1
                # bounded_positive_slope_slices.append(sl)
    return (n_regions,n_bounded,n_bounded_positive,
            n_unbounded,n_unbounded_positive)
            #bounded_positive_slope_slices,unbounded_positive_slope_slices)


def construct_poly(sname):
    return RNAPolytope.construct_from_file(sname)

CA1={}
CB1={}
CC1={}
for n in sobj_names:
    P = construct_poly(n)
    CA1[n],CB1[n],CC1[n]=get_bounded_region_intervals(P,nx1filter)

[2017-08-21 07:06:33.607218] [get_bounde] Start [2017-08-21 07:06:33.662287] [get_bounde] Stop [2017-08-21 07:06:46.740194] [get_bounde] Start [2017-08-21 07:06:46.820065] [get_bounde] Stop [2017-08-21 07:06:49.443636] [get_bounde] Start [2017-08-21 07:06:49.538779] [get_bounde] Stop [2017-08-21 07:06:52.485384] [get_bounde] Start [2017-08-21 07:06:52.600057] [get_bounde] Stop [2017-08-21 07:06:54.814969] [get_bounde] Start [2017-08-21 07:06:54.879715] [get_bounde] Stop [2017-08-21 07:06:56.981621] [get_bounde] Start [2017-08-21 07:06:57.046336] [get_bounde] Stop [2017-08-21 07:07:09.955820] [get_bounde] Start [2017-08-21 07:07:10.032743] [get_bounde] Stop [2017-08-21 07:07:12.970741] [get_bounde] Start [2017-08-21 07:07:13.053651] [get_bounde] Stop [2017-08-21 07:07:15.431441] [get_bounde] Start [2017-08-21 07:07:15.525095] [get_bounde] Stop [2017-08-21 07:07:18.422921] [get_bounde] Start [2017-08-21 07:07:18.524702] [get_bounde] Stop [2017-08-21 07:07:21.744981] [get_bounde] Start [2017-08-21 07:07:21.827369] [get_bounde] Stop [2017-08-21 07:07:35.060297] [get_bounde] Start [2017-08-21 07:07:35.175834] [get_bounde] Stop [2017-08-21 07:07:37.872059] [get_bounde] Start [2017-08-21 07:07:37.960732] [get_bounde] Stop [2017-08-21 07:07:40.228754] [get_bounde] Start [2017-08-21 07:07:40.310293] [get_bounde] Stop [2017-08-21 07:07:42.763657] [get_bounde] Start [2017-08-21 07:07:42.859867] [get_bounde] Stop [2017-08-21 07:07:45.175699] [get_bounde] Start [2017-08-21 07:07:45.246788] [get_bounde] Stop [2017-08-21 07:07:58.614223] [get_bounde] Start [2017-08-21 07:07:58.727234] [get_bounde] Stop [2017-08-21 07:08:00.844859] [get_bounde] Start [2017-08-21 07:08:00.912768] [get_bounde] Stop [2017-08-21 07:08:02.775255] [get_bounde] Start [2017-08-21 07:08:02.831655] [get_bounde] Stop [2017-08-21 07:08:05.316467] [get_bounde] Start [2017-08-21 07:08:05.410843] [get_bounde] Stop [2017-08-21 07:08:07.385700] [get_bounde] Start [2017-08-21 07:08:07.450675] [get_bounde] Stop [2017-08-21 07:08:20.711263] [get_bounde] Start [2017-08-21 07:08:20.802725] [get_bounde] Stop [2017-08-21 07:08:23.183542] [get_bounde] Start [2017-08-21 07:08:23.258708] [get_bounde] Stop [2017-08-21 07:08:25.679738] [get_bounde] Start [2017-08-21 07:08:25.760389] [get_bounde] Stop [2017-08-21 07:08:27.820120] [get_bounde] Start [2017-08-21 07:08:27.896698] [get_bounde] Stop [2017-08-21 07:08:30.697009] [get_bounde] Start [2017-08-21 07:08:30.775228] [get_bounde] Stop [2017-08-21 07:08:45.203333] [get_bounde] Start [2017-08-21 07:08:45.261480] [get_bounde] Stop [2017-08-21 07:08:48.411390] [get_bounde] Start [2017-08-21 07:08:48.510335] [get_bounde] Stop [2017-08-21 07:08:51.578302] [get_bounde] Start [2017-08-21 07:08:51.679712] [get_bounde] Stop [2017-08-21 07:08:53.916667] [get_bounde] Start [2017-08-21 07:08:54.005720] [get_bounde] Stop [2017-08-21 07:08:56.429874] [get_bounde] Start [2017-08-21 07:08:56.513454] [get_bounde] Stop [2017-08-21 07:09:10.838551] [get_bounde] Start [2017-08-21 07:09:10.918806] [get_bounde] Stop [2017-08-21 07:09:13.810443] [get_bounde] Start [2017-08-21 07:09:13.909707] [get_bounde] Stop [2017-08-21 07:09:16.394285] [get_bounde] Start [2017-08-21 07:09:16.510549] [get_bounde] Stop [2017-08-21 07:09:19.136141] [get_bounde] Start [2017-08-21 07:09:19.235452] [get_bounde] Stop [2017-08-21 07:09:21.275287] [get_bounde] Start [2017-08-21 07:09:21.346521] [get_bounde] Stop [2017-08-21 07:09:34.841336] [get_bounde] Start [2017-08-21 07:09:34.895589] [get_bounde] Stop [2017-08-21 07:09:38.197651] [get_bounde] Start [2017-08-21 07:09:38.304901] [get_bounde] Stop [2017-08-21 07:09:40.540849] [get_bounde] Start [2017-08-21 07:09:40.621749] [get_bounde] Stop

float(min([x[0] for n in CA1 for x in CA1[n]]))
float(max([x[1] for n in CA1 for x in CA1[n]]))
float(min([x[0] for n in CB1 for x in CB1[n]]))
float(max([x[1] for n in CB1 for x in CB1[n]]))
float(min([x[0] for n in CC1 for x in CC1[n]]))
float(max([x[1] for n in CC1 for x in CC1[n]]))

-1.55 86.8 -7.0 6.81 -52.7 -0.95

%md
# Box dimensions for bounded regions with x!= 1
```-1.55
86.8
-7.0
6.81
-52.7
-0.95```

Box dimensions for bounded regions with x!= 1

86.8
-7.0
6.81
-52.7
-0.95```

CA={}
CB={}
CC={}
for n in sobj_names:
    P = construct_poly(n)
    CA[n],CB[n],CC[n]=get_bounded_region_intervals(P)

[2017-08-21 06:19:21.151328] [get_bounde] Start [2017-08-21 06:19:21.243433] [get_bounde] Stop [2017-08-21 06:19:23.544423] [get_bounde] Start [2017-08-21 06:19:23.648140] [get_bounde] Stop [2017-08-21 06:19:26.295015] [get_bounde] Start [2017-08-21 06:19:26.386953] [get_bounde] Stop [2017-08-21 06:19:29.082228] [get_bounde] Start [2017-08-21 06:19:29.180020] [get_bounde] Stop [2017-08-21 06:19:31.390826] [get_bounde] Start [2017-08-21 06:19:31.478028] [get_bounde] Stop [2017-08-21 06:19:44.549958] [get_bounde] Start [2017-08-21 06:19:44.642948] [get_bounde] Stop [2017-08-21 06:19:46.946556] [get_bounde] Start [2017-08-21 06:19:47.043784] [get_bounde] Stop [2017-08-21 06:19:49.454277] [get_bounde] Start [2017-08-21 06:19:49.561823] [get_bounde] Stop [2017-08-21 06:19:51.839712] [get_bounde] Start [2017-08-21 06:19:51.926903] [get_bounde] Stop [2017-08-21 06:19:54.201189] [get_bounde] Start [2017-08-21 06:19:54.274265] [get_bounde] Stop [2017-08-21 06:20:07.273926] [get_bounde] Start [2017-08-21 06:20:07.347440] [get_bounde] Stop [2017-08-21 06:20:09.861869] [get_bounde] Start [2017-08-21 06:20:09.977792] [get_bounde] Stop [2017-08-21 06:20:12.268361] [get_bounde] Start [2017-08-21 06:20:12.365895] [get_bounde] Stop [2017-08-21 06:20:14.403863] [get_bounde] Start [2017-08-21 06:20:14.496030] [get_bounde] Stop [2017-08-21 06:20:16.678642] [get_bounde] Start [2017-08-21 06:20:16.765954] [get_bounde] Stop [2017-08-21 06:20:30.166802] [get_bounde] Start [2017-08-21 06:20:30.242432] [get_bounde] Stop [2017-08-21 06:20:33.209703] [get_bounde] Start [2017-08-21 06:20:33.332664] [get_bounde] Stop [2017-08-21 06:20:35.414301] [get_bounde] Start [2017-08-21 06:20:35.513596] [get_bounde] Stop [2017-08-21 06:20:37.308355] [get_bounde] Start [2017-08-21 06:20:37.378691] [get_bounde] Stop [2017-08-21 06:20:39.836041] [get_bounde] Start [2017-08-21 06:20:39.931071] [get_bounde] Stop [2017-08-21 06:20:41.902746] [get_bounde] Start [2017-08-21 06:20:41.979376] [get_bounde] Stop [2017-08-21 06:20:55.317550] [get_bounde] Start [2017-08-21 06:20:55.416869] [get_bounde] Stop [2017-08-21 06:20:57.593475] [get_bounde] Start [2017-08-21 06:20:57.681888] [get_bounde] Stop [2017-08-21 06:20:59.898745] [get_bounde] Start [2017-08-21 06:20:59.985440] [get_bounde] Stop [2017-08-21 06:21:01.983804] [get_bounde] Start [2017-08-21 06:21:02.068661] [get_bounde] Stop [2017-08-21 06:21:04.517776] [get_bounde] Start [2017-08-21 06:21:04.608932] [get_bounde] Stop [2017-08-21 06:21:18.325384] [get_bounde] Start [2017-08-21 06:21:18.435657] [get_bounde] Stop [2017-08-21 06:21:20.927682] [get_bounde] Start [2017-08-21 06:21:21.036601] [get_bounde] Stop [2017-08-21 06:21:23.415371] [get_bounde] Start [2017-08-21 06:21:23.516575] [get_bounde] Stop [2017-08-21 06:21:25.453224] [get_bounde] Start [2017-08-21 06:21:25.532705] [get_bounde] Stop [2017-08-21 06:21:27.829028] [get_bounde] Start [2017-08-21 06:21:27.923034] [get_bounde] Stop [2017-08-21 06:21:41.508231] [get_bounde] Start [2017-08-21 06:21:41.613213] [get_bounde] Stop [2017-08-21 06:21:44.111953] [get_bounde] Start [2017-08-21 06:21:44.223119] [get_bounde] Stop [2017-08-21 06:21:46.404373] [get_bounde] Start [2017-08-21 06:21:46.497411] [get_bounde] Stop [2017-08-21 06:21:48.854676] [get_bounde] Start [2017-08-21 06:21:48.940207] [get_bounde] Stop [2017-08-21 06:21:50.926142] [get_bounde] Start [2017-08-21 06:21:51.009006] [get_bounde] Stop [2017-08-21 06:22:03.968055] [get_bounde] Start [2017-08-21 06:22:04.066633] [get_bounde] Stop [2017-08-21 06:22:06.669229] [get_bounde] Start [2017-08-21 06:22:06.778215] [get_bounde] Stop [2017-08-21 06:22:08.653081] [get_bounde] Start [2017-08-21 06:22:08.741599] [get_bounde] Stop

float(min([x[0] for n in CA for x in CA[n]]))
float(max([x[1] for n in CA for x in CA[n]]))
float(min([x[0] for n in CB for x in CB[n]]))
float(max([x[1] for n in CB for x in CB[n]]))
float(min([x[0] for n in CC for x in CC[n]]))
float(max([x[1] for n in CC for x in CC[n]]))

-11.478947368421053 540.75 -7.0 7.1 -52.7 5.089473684210526

%md
# Box dimensions for bounded regions with x= 1
```
-11.478947368421053
540.75
-7.0
7.1
-52.7
5.089473684210526
```

Box dimensions for bounded regions with x= 1

-11.478947368421053
540.75
-7.0
7.1
-52.7
5.089473684210526

# 3 mins
A1m={}
B1m={}
C1m={}
x1filter = lambda sl: sl.original_vertex.vector()[0]>1
for n in sobj_names:
    P=construct_poly(n)
    A1m[n],B1m[n],C1m[n] = get_bounded_region_widths(P,x1filter)

[2017-08-21 05:33:51.200582] [get_bounde] Start [2017-08-21 05:33:51.268148] [get_bounde] Stop [2017-08-21 05:33:53.318064] [get_bounde] Start [2017-08-21 05:33:53.378157] [get_bounde] Stop [2017-08-21 05:34:06.516670] [get_bounde] Start [2017-08-21 05:34:06.573566] [get_bounde] Stop [2017-08-21 05:34:09.366793] [get_bounde] Start [2017-08-21 05:34:09.477880] [get_bounde] Stop [2017-08-21 05:34:12.010714] [get_bounde] Start [2017-08-21 05:34:12.085254] [get_bounde] Stop [2017-08-21 05:34:14.238231] [get_bounde] Start [2017-08-21 05:34:14.307175] [get_bounde] Stop [2017-08-21 05:34:16.349958] [get_bounde] Start [2017-08-21 05:34:16.413936] [get_bounde] Stop [2017-08-21 05:34:29.240635] [get_bounde] Start [2017-08-21 05:34:29.287641] [get_bounde] Stop [2017-08-21 05:34:31.865734] [get_bounde] Start [2017-08-21 05:34:31.957799] [get_bounde] Stop [2017-08-21 05:34:34.471270] [get_bounde] Start [2017-08-21 05:34:34.558712] [get_bounde] Stop [2017-08-21 05:34:36.900423] [get_bounde] Start [2017-08-21 05:34:36.977360] [get_bounde] Stop [2017-08-21 05:34:39.277713] [get_bounde] Start [2017-08-21 05:34:39.343002] [get_bounde] Stop [2017-08-21 05:34:53.140607] [get_bounde] Start [2017-08-21 05:34:53.184327] [get_bounde] Stop [2017-08-21 05:34:55.354042] [get_bounde] Start [2017-08-21 05:34:55.443568] [get_bounde] Stop [2017-08-21 05:34:57.919495] [get_bounde] Start [2017-08-21 05:34:58.014978] [get_bounde] Stop [2017-08-21 05:35:00.469221] [get_bounde] Start [2017-08-21 05:35:00.551452] [get_bounde] Stop [2017-08-21 05:35:03.286886] [get_bounde] Start [2017-08-21 05:35:03.360541] [get_bounde] Stop [2017-08-21 05:35:05.241227] [get_bounde] Start [2017-08-21 05:35:05.303873] [get_bounde] Stop [2017-08-21 05:35:17.607569] [get_bounde] Start [2017-08-21 05:35:17.661298] [get_bounde] Stop [2017-08-21 05:35:20.129088] [get_bounde] Start [2017-08-21 05:35:20.224464] [get_bounde] Stop [2017-08-21 05:35:22.186109] [get_bounde] Start [2017-08-21 05:35:22.252276] [get_bounde] Stop [2017-08-21 05:35:24.468983] [get_bounde] Start [2017-08-21 05:35:24.545471] [get_bounde] Stop [2017-08-21 05:35:26.559309] [get_bounde] Start [2017-08-21 05:35:26.625315] [get_bounde] Stop [2017-08-21 05:35:39.194392] [get_bounde] Start [2017-08-21 05:35:39.250263] [get_bounde] Stop [2017-08-21 05:35:41.305579] [get_bounde] Start [2017-08-21 05:35:41.378342] [get_bounde] Stop [2017-08-21 05:35:43.835112] [get_bounde] Start [2017-08-21 05:35:43.916130] [get_bounde] Stop [2017-08-21 05:35:46.362811] [get_bounde] Start [2017-08-21 05:35:46.445660] [get_bounde] Stop [2017-08-21 05:35:48.697667] [get_bounde] Start [2017-08-21 05:35:48.774300] [get_bounde] Stop [2017-08-21 05:36:01.613755] [get_bounde] Start [2017-08-21 05:36:01.679035] [get_bounde] Stop [2017-08-21 05:36:03.712812] [get_bounde] Start [2017-08-21 05:36:03.782442] [get_bounde] Stop [2017-08-21 05:36:06.134082] [get_bounde] Start [2017-08-21 05:36:06.223886] [get_bounde] Stop [2017-08-21 05:36:08.649394] [get_bounde] Start [2017-08-21 05:36:08.726165] [get_bounde] Stop [2017-08-21 05:36:11.085094] [get_bounde] Start [2017-08-21 05:36:11.156734] [get_bounde] Stop [2017-08-21 05:36:23.871203] [get_bounde] Start [2017-08-21 05:36:23.965252] [get_bounde] Stop [2017-08-21 05:36:26.358499] [get_bounde] Start [2017-08-21 05:36:26.453427] [get_bounde] Stop [2017-08-21 05:36:28.346117] [get_bounde] Start [2017-08-21 05:36:28.419940] [get_bounde] Stop [2017-08-21 05:36:30.520314] [get_bounde] Start [2017-08-21 05:36:30.585950] [get_bounde] Stop [2017-08-21 05:36:32.987152] [get_bounde] Start [2017-08-21 05:36:33.065532] [get_bounde] Stop [2017-08-21 05:36:45.241474] [get_bounde] Start [2017-08-21 05:36:45.302834] [get_bounde] Stop

A1ms = []
B1ms = []
C1ms = []
for n in A:
    A1ms.extend(A1m[n])
    B1ms.extend(B1m[n])
    C1ms.extend(C1m[n])
avg(A1ms)
std(A1ms)
avg(B1ms)
std(B1ms)
avg(C1ms)
std(C1ms)

5.1547776862131949 4.8535618051007487 0.45730904256640104 0.50965601903647251 2.6985742681351845 2.7747218357039536

float(max(A1ms))
float(min(A1ms))
float(max(B1ms))
float(min(B1ms))
float(max(C1ms))
float(min(C1ms))

51.2 0.03 6.076666666666667 0.0025 30.9 0.018333333333333333

# 3 mins
A1={}
B1={}
C1={}
x1filter = lambda sl: sl.original_vertex.vector()[0]==1
for n in sobj_names:
    P=construct_poly(n)
    A1[n],B1[n],C1[n] = get_bounded_region_widths(P,x1filter)

[2017-08-21 05:24:21.974304] [get_bounde] Start [2017-08-21 05:24:21.990461] [get_bounde] Stop [2017-08-21 05:24:35.661086] [get_bounde] Start [2017-08-21 05:24:35.704295] [get_bounde] Stop [2017-08-21 05:24:38.234997] [get_bounde] Start [2017-08-21 05:24:38.258556] [get_bounde] Stop [2017-08-21 05:24:40.960114] [get_bounde] Start [2017-08-21 05:24:40.979737] [get_bounde] Stop [2017-08-21 05:24:43.243925] [get_bounde] Start [2017-08-21 05:24:43.266264] [get_bounde] Stop [2017-08-21 05:24:45.387771] [get_bounde] Start [2017-08-21 05:24:45.417618] [get_bounde] Stop [2017-08-21 05:24:58.456715] [get_bounde] Start [2017-08-21 05:24:58.478252] [get_bounde] Stop [2017-08-21 05:25:00.964309] [get_bounde] Start [2017-08-21 05:25:01.002777] [get_bounde] Stop [2017-08-21 05:25:03.266503] [get_bounde] Start [2017-08-21 05:25:03.295711] [get_bounde] Stop [2017-08-21 05:25:05.573751] [get_bounde] Start [2017-08-21 05:25:05.602546] [get_bounde] Stop [2017-08-21 05:25:07.957757] [get_bounde] Start [2017-08-21 05:25:07.976569] [get_bounde] Stop [2017-08-21 05:25:21.033889] [get_bounde] Start [2017-08-21 05:25:21.063890] [get_bounde] Stop [2017-08-21 05:25:23.298191] [get_bounde] Start [2017-08-21 05:25:23.326958] [get_bounde] Stop [2017-08-21 05:25:25.343198] [get_bounde] Start [2017-08-21 05:25:25.362933] [get_bounde] Stop [2017-08-21 05:25:27.602004] [get_bounde] Start [2017-08-21 05:25:27.621645] [get_bounde] Stop [2017-08-21 05:25:30.169497] [get_bounde] Start [2017-08-21 05:25:30.202102] [get_bounde] Stop [2017-08-21 05:25:43.695554] [get_bounde] Start [2017-08-21 05:25:43.716633] [get_bounde] Stop [2017-08-21 05:25:45.667357] [get_bounde] Start [2017-08-21 05:25:45.689178] [get_bounde] Stop [2017-08-21 05:25:47.411500] [get_bounde] Start [2017-08-21 05:25:47.431929] [get_bounde] Stop [2017-08-21 05:25:49.813893] [get_bounde] Start [2017-08-21 05:25:49.833578] [get_bounde] Stop [2017-08-21 05:25:51.945314] [get_bounde] Start [2017-08-21 05:25:51.964700] [get_bounde] Stop [2017-08-21 05:26:05.016153] [get_bounde] Start [2017-08-21 05:26:05.044208] [get_bounde] Stop [2017-08-21 05:26:07.507487] [get_bounde] Start [2017-08-21 05:26:07.532535] [get_bounde] Stop [2017-08-21 05:26:10.253984] [get_bounde] Start [2017-08-21 05:26:10.273763] [get_bounde] Stop [2017-08-21 05:26:12.308449] [get_bounde] Start [2017-08-21 05:26:12.327467] [get_bounde] Stop [2017-08-21 05:26:15.343686] [get_bounde] Start [2017-08-21 05:26:15.364920] [get_bounde] Stop [2017-08-21 05:26:28.799085] [get_bounde] Start [2017-08-21 05:26:28.824919] [get_bounde] Stop [2017-08-21 05:26:31.306505] [get_bounde] Start [2017-08-21 05:26:31.332263] [get_bounde] Stop [2017-08-21 05:26:33.726147] [get_bounde] Start [2017-08-21 05:26:33.748548] [get_bounde] Stop [2017-08-21 05:26:35.716238] [get_bounde] Start [2017-08-21 05:26:35.741717] [get_bounde] Stop [2017-08-21 05:26:38.146236] [get_bounde] Start [2017-08-21 05:26:38.169716] [get_bounde] Stop [2017-08-21 05:26:51.460247] [get_bounde] Start [2017-08-21 05:26:51.485865] [get_bounde] Stop [2017-08-21 05:26:54.117507] [get_bounde] Start [2017-08-21 05:26:54.148034] [get_bounde] Stop [2017-08-21 05:26:56.547467] [get_bounde] Start [2017-08-21 05:26:56.568652] [get_bounde] Stop [2017-08-21 05:26:59.132729] [get_bounde] Start [2017-08-21 05:26:59.156911] [get_bounde] Stop [2017-08-21 05:27:01.034524] [get_bounde] Start [2017-08-21 05:27:14.834315] [get_bounde] Stop [2017-08-21 05:27:17.584735] [get_bounde] Start [2017-08-21 05:27:17.612033] [get_bounde] Stop [2017-08-21 05:27:19.675486] [get_bounde] Start [2017-08-21 05:27:19.702898] [get_bounde] Stop

A1s = []
B1s = []
C1s = []
for n in A:
    A1s.extend(A1[n])
    B1s.extend(B1[n])
    C1s.extend(C1[n])
avg(A1s)
std(A1s)
avg(B1s)
std(B1s)
avg(C1s)
std(C1s)

84.938675331207975 64.864962492033541 0.39990983036369532 0.51242429289672409 2.3815814021630222 2.3901260145979997

len(A1s)

3608

# 3 mins
A={}
B={}
C={}
for n in sobj_names:
    P=construct_poly(n)
    A[n],B[n],C[n] = get_bounded_region_widths(P)

[2017-08-21 04:09:14.507355] [get_bounde] Start [2017-08-21 04:09:14.572727] [get_bounde] Stop [2017-08-21 04:09:16.611765] [get_bounde] Start [2017-08-21 04:09:16.677078] [get_bounde] Stop [2017-08-21 04:09:19.003670] [get_bounde] Start [2017-08-21 04:09:19.081420] [get_bounde] Stop [2017-08-21 04:09:21.662772] [get_bounde] Start [2017-08-21 04:09:21.757646] [get_bounde] Stop [2017-08-21 04:09:35.923343] [get_bounde] Start [2017-08-21 04:09:36.020389] [get_bounde] Stop [2017-08-21 04:09:38.042732] [get_bounde] Start [2017-08-21 04:09:38.129541] [get_bounde] Stop [2017-08-21 04:09:40.162388] [get_bounde] Start [2017-08-21 04:09:40.241690] [get_bounde] Stop [2017-08-21 04:09:42.640886] [get_bounde] Start [2017-08-21 04:09:42.727342] [get_bounde] Stop [2017-08-21 04:09:45.083307] [get_bounde] Start [2017-08-21 04:09:45.159156] [get_bounde] Stop [2017-08-21 04:09:58.546610] [get_bounde] Start [2017-08-21 04:09:58.622945] [get_bounde] Stop [2017-08-21 04:10:01.040043] [get_bounde] Start [2017-08-21 04:10:01.146141] [get_bounde] Stop [2017-08-21 04:10:03.768182] [get_bounde] Start [2017-08-21 04:10:03.861339] [get_bounde] Stop [2017-08-21 04:10:06.207852] [get_bounde] Start [2017-08-21 04:10:06.277789] [get_bounde] Stop [2017-08-21 04:10:08.341236] [get_bounde] Start [2017-08-21 04:10:08.414089] [get_bounde] Stop [2017-08-21 04:10:21.321264] [get_bounde] Start [2017-08-21 04:10:21.370515] [get_bounde] Stop [2017-08-21 04:10:23.899122] [get_bounde] Start [2017-08-21 04:10:24.006065] [get_bounde] Stop [2017-08-21 04:10:26.865495] [get_bounde] Start [2017-08-21 04:10:26.986512] [get_bounde] Stop [2017-08-21 04:10:29.068746] [get_bounde] Start [2017-08-21 04:10:29.148663] [get_bounde] Stop [2017-08-21 04:10:30.941192] [get_bounde] Start [2017-08-21 04:10:31.013449] [get_bounde] Stop [2017-08-21 04:10:33.524748] [get_bounde] Start [2017-08-21 04:10:33.616210] [get_bounde] Stop [2017-08-21 04:10:46.244616] [get_bounde] Start [2017-08-21 04:10:46.348423] [get_bounde] Stop [2017-08-21 04:10:48.680293] [get_bounde] Start [2017-08-21 04:10:48.769631] [get_bounde] Stop [2017-08-21 04:10:50.824248] [get_bounde] Start [2017-08-21 04:10:50.896799] [get_bounde] Stop [2017-08-21 04:10:53.062181] [get_bounde] Start [2017-08-21 04:10:53.148882] [get_bounde] Stop [2017-08-21 04:10:54.983642] [get_bounde] Start [2017-08-21 04:10:55.044047] [get_bounde] Stop [2017-08-21 04:11:07.717962] [get_bounde] Start [2017-08-21 04:11:07.827819] [get_bounde] Stop [2017-08-21 04:11:10.537733] [get_bounde] Start [2017-08-21 04:11:10.643549] [get_bounde] Stop [2017-08-21 04:11:13.169401] [get_bounde] Start [2017-08-21 04:11:13.267798] [get_bounde] Stop [2017-08-21 04:11:15.614652] [get_bounde] Start [2017-08-21 04:11:15.700996] [get_bounde] Stop [2017-08-21 04:11:17.701565] [get_bounde] Start [2017-08-21 04:11:17.770926] [get_bounde] Stop [2017-08-21 04:11:31.177551] [get_bounde] Start [2017-08-21 04:11:31.274111] [get_bounde] Stop [2017-08-21 04:11:33.942489] [get_bounde] Start [2017-08-21 04:11:34.045612] [get_bounde] Stop [2017-08-21 04:11:36.505916] [get_bounde] Start [2017-08-21 04:11:36.599725] [get_bounde] Stop [2017-08-21 04:11:38.871916] [get_bounde] Start [2017-08-21 04:11:38.968902] [get_bounde] Stop [2017-08-21 04:11:41.297335] [get_bounde] Start [2017-08-21 04:11:41.379737] [get_bounde] Stop [2017-08-21 04:11:53.801067] [get_bounde] Start [2017-08-21 04:11:53.880181] [get_bounde] Stop [2017-08-21 04:11:56.099921] [get_bounde] Start [2017-08-21 04:11:56.184983] [get_bounde] Stop [2017-08-21 04:11:58.779081] [get_bounde] Start [2017-08-21 04:11:58.875520] [get_bounde] Stop [2017-08-21 04:12:00.689100] [get_bounde] Start [2017-08-21 04:12:00.761975] [get_bounde] Stop

As = []
Bs = []
Cs = []
for n in A:
    As.extend(A[n])
    Bs.extend(B[n])
    Cs.extend(C[n])
avg(As)
std(As)
avg(Bs)
std(Bs)
avg(Cs)
std(Cs)

28.259385014501273 50.445283479276704 0.44068681304337654 0.51112271239128726 2.6067763501299224 2.6729232154816196

positive_orthant_sigs = {}
volumes = []
print timestamp()
for sname in sobj_names:
    print timestamp(sname)
    P = construct_poly(sname)
    vol,positive_orthant_sigs[sname] = collect_positive_orthant_signatures(P)
    volumes.append(vol)
print timestamp()

[2017-08-21 00:23:14.044761] [2017-08-21 00:23:14.045326] sobjs/nitrosomonas.europaea.atcc.19718_BX321863_b.sobj [2017-08-21 00:23:18.667964] sobjs/lactobacillus.casei.atcc.334_CP000423_j.sobj [2017-08-21 00:23:24.184594] sobjs/methanosarcina.barkeri.str..fusaro_CP000099_bi.sobj [2017-08-21 00:23:28.783189] sobjs/uncultured.alpha.proteobacterium_AY933788.sobj [2017-08-21 00:23:35.673521] sobjs/tropheryma.whipplei.str..twist_AE014184_h.sobj [2017-08-21 00:23:41.389837] sobjs/uncultured.thiotrichales.bacterium_AY934254_a.sobj [2017-08-21 00:23:46.496196] sobjs/synechococcus.sp..wh.8102_BX569691.sobj [2017-08-21 00:23:52.952927] sobjs/yersinia.pestis.kim_AE009952_ae.sobj [2017-08-21 00:23:58.474454] sobjs/uncultured.alpha.proteobacterium_DQ093144_a.sobj [2017-08-21 00:24:04.545777] sobjs/carboxydothermus.hydrogenoformans.z-2901_CP000141_af.sobj [2017-08-21 00:24:11.319464] sobjs/haloferax.volcanii_HALTRK1.sobj [2017-08-21 00:24:14.906140] sobjs/homo.sapiens_X04779.sobj [2017-08-21 00:24:22.024994] sobjs/arabidopsis.thaliana_AC006340.sobj [2017-08-21 00:24:28.982465] sobjs/chlamydophila.pneumoniae.tw-183_AE017159_j.sobj [2017-08-21 00:24:32.529490] sobjs/thermoplasma.volcanium.gss1_BA000011_b.sobj [2017-08-21 00:24:39.288231] sobjs/uncultured.gamma.proteobacterium_AY934184_a.sobj [2017-08-21 00:24:43.237021] sobjs/nocardia.farcinica.ifm.10152_AP006618_az.sobj [2017-08-21 00:24:51.238808] sobjs/chlorobium.phaeobacteroides.dsm.266_CP000492_ab.sobj [2017-08-21 00:24:54.900037] sobjs/scenedesmus.obliquus_DQ396875_n.sobj [2017-08-21 00:25:04.307624] sobjs/pyrococcus.abyssi.ge5_AJ248288.sobj [2017-08-21 00:25:09.748096] sobjs/euglena.longa_AJ294725_s.sobj [2017-08-21 00:25:13.865912] sobjs/marchantia.polymorpha_M.CAU.C.Mar.polym.1.sobj [2017-08-21 00:25:22.735115] sobjs/lactobacillus.delbrueckii.subsp..bulgaricus.atcc.baa-365_CP000412_cj.sobj [2017-08-21 00:25:28.280689] sobjs/uncultured.alteromonadales.bacterium_AY934393.sobj [2017-08-21 00:25:32.354818] sobjs/methanospirillum.hungatei.jf-1_CP000254_ad.sobj [2017-08-21 00:25:41.396964] sobjs/drosophila.melanogaster_DROTRF2.sobj [2017-08-21 00:25:46.119669] sobjs/pelobacter.carbinolicus.dsm.2380_CP000142_af.sobj [2017-08-21 00:25:51.132704] sobjs/bordetella.parapertussis_BX640433_e.sobj [2017-08-21 00:26:03.244447] sobjs/magnetococcus.sp..mc-1_CP000471_ao.sobj [2017-08-21 00:26:07.878783] sobjs/yarrowia.lipolytica.clib122_CR382129_by.sobj [2017-08-21 00:26:12.600207] sobjs/hyperthermus.butylicus.dsm.5456_CP000493_ad.sobj [2017-08-21 00:26:24.667616] sobjs/homo.sapiens_AC004932_k.sobj [2017-08-21 00:26:29.347073] sobjs/uncultured.gamma.proteobacterium_AY934387_a.sobj [2017-08-21 00:26:34.678523] sobjs/rickettsia.conorii.str..malish.7_AE008623.sobj [2017-08-21 00:26:38.853285] sobjs/solibacter.usitatus.ellin6076_CP000473_ay.sobj [2017-08-21 00:26:51.441757] sobjs/wigglesworthia.glossinidia.endosymbiont.of.glossina.brevipalpis_BA000021_ad.sobj [2017-08-21 00:26:56.082593] sobjs/sulfolobus.tokodaii.str..7_BA000023_ak.sobj [2017-08-21 00:27:00.394430] sobjs/uncultured.marinobacter.sp._AY934351.sobj [2017-08-21 00:27:05.019680] sobjs/nephroselmis.olivacea_AF137379_g.sobj [2017-08-21 00:27:18.470749]

avg(volumes)
std(volumes)

3.2556192400903194 2.8228455358801376

one_four_counter =0
for n in positive_orthant_sigs:
    print "{0}: {1}".format(n,positive_orthant_sigs[n].keys())
    if (1,4) in positive_orthant_sigs[n]:
        one_four_counter+=1
    for k in positive_orthant_sigs[n]:
        for sl in positive_orthant_sigs[n][k]:
            x,y,z,w = sl.original_vertex.vector()
            print "{0}, ".format((x,y,z)),
        print ""
print one_four_counter
print len(positive_orthant_sigs.keys())

sobjs/methanosarcina.barkeri.str..fusaro_CP000099_bi.sobj: [(0, 0), (1, 4)] (0, 0, 0), (1, 6, 4), sobjs/chlorobium.phaeobacteroides.dsm.266_CP000492_ab.sobj: [(0, 0), (1, 3), (1, 4)] (0, 0, 0), (1, 3, 3), (1, 1, 3), (1, 13, 4), sobjs/homo.sapiens_X04779.sobj: [(0, 0), (2, 6), (1, 3), (1, 4)] (0, 0, 0), (2, 7, 6), (2, 6, 6), (1, 2, 3), (1, 5, 3), (1, 8, 4), sobjs/solibacter.usitatus.ellin6076_CP000473_ay.sobj: [(1, 3), (0, 0), (2, 6), (1, 4)] (1, 8, 3), (0, 0, 0), (2, 9, 6), (2, 4, 6), (1, 8, 4), sobjs/uncultured.gamma.proteobacterium_AY934387_a.sobj: [(1, 3), (2, 6), (0, 0)] (1, 2, 3), (1, 7, 3), (1, 11, 3), (2, 7, 6), (2, 10, 6), (0, 0, 0), sobjs/pyrococcus.abyssi.ge5_AJ248288.sobj: [(1, 3), (2, 6), (0, 0), (1, 4)] (1, 2, 3), (1, 9, 3), (2, 4, 6), (2, 1, 6), (2, 6, 6), (2, 8, 6), (0, 0, 0), (1, 6, 4), sobjs/yersinia.pestis.kim_AE009952_ae.sobj: [(1, 3), (0, 0)] (1, 3, 3), (1, 5, 3), (0, 0, 0), sobjs/yarrowia.lipolytica.clib122_CR382129_by.sobj: [(2, 6), (1, 3), (0, 0)] (2, 7, 6), (2, 5, 6), (2, 4, 6), (1, 1, 3), (1, 10, 3), (1, 6, 3), (0, 0, 0), sobjs/sulfolobus.tokodaii.str..7_BA000023_ak.sobj: [(0, 0), (1, 3), (1, 4)] (0, 0, 0), (1, 18, 3), (1, 3, 4), (1, 1, 4), sobjs/homo.sapiens_AC004932_k.sobj: [(1, 3), (2, 6), (0, 0)] (1, 2, 3), (1, 7, 3), (2, 7, 6), (2, 12, 6), (0, 0, 0), sobjs/euglena.longa_AJ294725_s.sobj: [(1, 3), (0, 0), (2, 6), (1, 4)] (1, 0, 3), (1, 22, 3), (0, 0, 0), (2, 4, 6), (1, 8, 4), sobjs/bordetella.parapertussis_BX640433_e.sobj: [(1, 3), (0, 0), (2, 6), (1, 4)] (1, 17, 3), (0, 0, 0), (2, 3, 6), (2, 5, 6), (2, 8, 6), (2, 10, 6), (1, 9, 4), (1, 1, 4), (1, 3, 4), sobjs/nitrosomonas.europaea.atcc.19718_BX321863_b.sobj: [(1, 3), (0, 0)] (1, 9, 3), (1, 5, 3), (1, 3, 3), (0, 0, 0), sobjs/rickettsia.conorii.str..malish.7_AE008623.sobj: [(2, 6), (0, 0), (1, 4)] (2, 6, 6), (0, 0, 0), (1, 0, 4), (1, 8, 4), sobjs/uncultured.alpha.proteobacterium_AY933788.sobj: [(1, 3), (2, 6), (0, 0), (1, 4)] (1, 3, 3), (2, 7, 6), (2, 9, 6), (2, 5, 6), (0, 0, 0), (1, 0, 4), (1, 8, 4), sobjs/methanospirillum.hungatei.jf-1_CP000254_ad.sobj: [(0, 0), (1, 3), (1, 4)] (0, 0, 0), (1, 19, 3), (1, 2, 4), (1, 4, 4), sobjs/wigglesworthia.glossinidia.endosymbiont.of.glossina.brevipalpis_BA000021_ad.sobj: [(1, 3), (2, 6), (0, 0), (1, 4)] (1, 17, 3), (2, 5, 6), (2, 7, 6), (0, 0, 0), (1, 2, 4), (1, 0, 4), sobjs/arabidopsis.thaliana_AC006340.sobj: [(1, 3), (2, 6), (0, 0), (1, 4)] (1, 18, 3), (2, 6, 6), (0, 0, 0), (1, 4, 4), (1, 9, 4), sobjs/synechococcus.sp..wh.8102_BX569691.sobj: [(1, 3), (2, 6), (1, 4), (0, 0)] (1, 1, 3), (1, 3, 3), (2, 2, 6), (2, 8, 6), (2, 6, 6), (2, 4, 6), (1, 6, 4), (0, 0, 0), sobjs/drosophila.melanogaster_DROTRF2.sobj: [(1, 3), (2, 6), (0, 0), (1, 4)] (1, 14, 3), (1, 8, 3), (2, 3, 6), (0, 0, 0), (1, 5, 4), sobjs/magnetococcus.sp..mc-1_CP000471_ao.sobj: [(1, 3), (2, 6), (0, 0), (1, 4)] (1, 3, 3), (2, 4, 6), (2, 10, 6), (0, 0, 0), (1, 6, 4), sobjs/pelobacter.carbinolicus.dsm.2380_CP000142_af.sobj: [(0, 0), (2, 6), (1, 3), (1, 4)] (0, 0, 0), (2, 6, 6), (1, 4, 3), (1, 6, 3), (1, 8, 4), sobjs/nephroselmis.olivacea_AF137379_g.sobj: [(2, 6), (1, 3), (1, 4), (0, 0)] (2, 7, 6), (2, 9, 6), (1, 19, 3), (1, 1, 3), (1, 3, 3), (1, 13, 4), (0, 0, 0), sobjs/tropheryma.whipplei.str..twist_AE014184_h.sobj: [(0, 0), (1, 4)] (0, 0, 0), (1, 2, 4), (1, 4, 4), sobjs/carboxydothermus.hydrogenoformans.z-2901_CP000141_af.sobj: [(0, 0), (2, 6), (1, 4)] (0, 0, 0), (2, 4, 6), (2, 6, 6), (1, 7, 4), (1, 2, 4), (1, 9, 4), sobjs/scenedesmus.obliquus_DQ396875_n.sobj: [(2, 6), (1, 3), (1, 4), (0, 0)] (2, 12, 6), (2, 7, 6), (2, 2, 6), (2, 9, 6), (1, 1, 3), (1, 4, 3), (1, 9, 3), (1, 6, 3), (1, 9, 4), (1, 11, 4), (0, 0, 0), sobjs/uncultured.gamma.proteobacterium_AY934184_a.sobj: [(1, 3), (0, 0)] (1, 3, 3), (1, 4, 3), (1, 8, 3), (0, 0, 0), sobjs/lactobacillus.casei.atcc.334_CP000423_j.sobj: [(2, 6), (1, 3), (0, 0)] (2, 16, 6), (1, 7, 3), (1, 5, 3), (1, 9, 3), (1, 3, 3), (0, 0, 0), sobjs/uncultured.alteromonadales.bacterium_AY934393.sobj: [(0, 0), (2, 6), (1, 4)] (0, 0, 0), (2, 5, 6), (1, 8, 4), (1, 0, 4), sobjs/uncultured.thiotrichales.bacterium_AY934254_a.sobj: [(1, 3), (0, 0)] (1, 7, 3), (1, 2, 3), (0, 0, 0), sobjs/uncultured.marinobacter.sp._AY934351.sobj: [(2, 6), (0, 0)] (2, 7, 6), (2, 9, 6), (2, 3, 6), (0, 0, 0), sobjs/uncultured.alpha.proteobacterium_DQ093144_a.sobj: [(1, 3), (0, 0), (2, 6), (1, 4)] (1, 3, 3), (1, 1, 3), (0, 0, 0), (2, 9, 6), (1, 6, 4), sobjs/chlamydophila.pneumoniae.tw-183_AE017159_j.sobj: [(2, 6), (0, 0), (1, 4)] (2, 7, 6), (2, 5, 6), (0, 0, 0), (1, 7, 4), sobjs/haloferax.volcanii_HALTRK1.sobj: [(0, 0), (1, 4)] (0, 0, 0), (1, 5, 4), (1, 3, 4), sobjs/thermoplasma.volcanium.gss1_BA000011_b.sobj: [(1, 3), (0, 0)] (1, 9, 3), (1, 6, 3), (0, 0, 0), sobjs/marchantia.polymorpha_M.CAU.C.Mar.polym.1.sobj: [(1, 3), (2, 6), (0, 0)] (1, 5, 3), (1, 3, 3), (2, 11, 6), (0, 0, 0), sobjs/nocardia.farcinica.ifm.10152_AP006618_az.sobj: [(2, 6), (1, 4), (0, 0)] (2, 1, 6), (1, 3, 4), (0, 0, 0), sobjs/hyperthermus.butylicus.dsm.5456_CP000493_ad.sobj: [(0, 0), (2, 6), (1, 4)] (0, 0, 0), (2, 13, 6), (2, 7, 6), (2, 6, 6), (1, 3, 4), (1, 9, 4), (1, 4, 4), sobjs/lactobacillus.delbrueckii.subsp..bulgaricus.atcc.baa-365_CP000412_cj.sobj: [(1, 3), (0, 0), (2, 6), (1, 4)] (1, 3, 3), (0, 0, 0), (2, 5, 6), (2, 8, 6), (2, 6, 6), (2, 10, 6), (2, 3, 6), (1, 8, 4), (1, 6, 4), 28 39

# ~4min
positive_orthant_slices = {}
print timestamp()
for sname in sobj_names:
    P = construct_poly(sname)
    positive_orthant_slices[sname] = count_positive_orthant_slices(P)
print timestamp()

for n in positive_orthant_slices:
    print "{0}: {1}".format(positive_orthant_slices[n],n)

3: sobjs/methanosarcina.barkeri.str..fusaro_CP000099_bi.sobj 5: sobjs/chlorobium.phaeobacteroides.dsm.266_CP000492_ab.sobj 6: sobjs/homo.sapiens_X04779.sobj 5: sobjs/solibacter.usitatus.ellin6076_CP000473_ay.sobj 6: sobjs/uncultured.gamma.proteobacterium_AY934387_a.sobj 9: sobjs/pyrococcus.abyssi.ge5_AJ248288.sobj 4: sobjs/yersinia.pestis.kim_AE009952_ae.sobj 10: sobjs/yarrowia.lipolytica.clib122_CR382129_by.sobj 4: sobjs/sulfolobus.tokodaii.str..7_BA000023_ak.sobj 5: sobjs/homo.sapiens_AC004932_k.sobj 5: sobjs/euglena.longa_AJ294725_s.sobj 9: sobjs/bordetella.parapertussis_BX640433_e.sobj 6: sobjs/nitrosomonas.europaea.atcc.19718_BX321863_b.sobj 4: sobjs/rickettsia.conorii.str..malish.7_AE008623.sobj 7: sobjs/uncultured.alpha.proteobacterium_AY933788.sobj 4: sobjs/methanospirillum.hungatei.jf-1_CP000254_ad.sobj 6: sobjs/wigglesworthia.glossinidia.endosymbiont.of.glossina.brevipalpis_BA000021_ad.sobj 9: sobjs/arabidopsis.thaliana_AC006340.sobj 8: sobjs/synechococcus.sp..wh.8102_BX569691.sobj 6: sobjs/drosophila.melanogaster_DROTRF2.sobj 6: sobjs/magnetococcus.sp..mc-1_CP000471_ao.sobj 5: sobjs/pelobacter.carbinolicus.dsm.2380_CP000142_af.sobj 7: sobjs/nephroselmis.olivacea_AF137379_g.sobj 3: sobjs/tropheryma.whipplei.str..twist_AE014184_h.sobj 8: sobjs/carboxydothermus.hydrogenoformans.z-2901_CP000141_af.sobj 12: sobjs/scenedesmus.obliquus_DQ396875_n.sobj 4: sobjs/uncultured.gamma.proteobacterium_AY934184_a.sobj 7: sobjs/lactobacillus.casei.atcc.334_CP000423_j.sobj 4: sobjs/uncultured.alteromonadales.bacterium_AY934393.sobj 4: sobjs/uncultured.thiotrichales.bacterium_AY934254_a.sobj 4: sobjs/uncultured.marinobacter.sp._AY934351.sobj 5: sobjs/uncultured.alpha.proteobacterium_DQ093144_a.sobj 4: sobjs/chlamydophila.pneumoniae.tw-183_AE017159_j.sobj 3: sobjs/haloferax.volcanii_HALTRK1.sobj 3: sobjs/thermoplasma.volcanium.gss1_BA000011_b.sobj 5: sobjs/marchantia.polymorpha_M.CAU.C.Mar.polym.1.sobj 4: sobjs/nocardia.farcinica.ifm.10152_AP006618_az.sobj 8: sobjs/hyperthermus.butylicus.dsm.5456_CP000493_ad.sobj 10: sobjs/lactobacillus.delbrueckii.subsp..bulgaricus.atcc.baa-365_CP000412_cj.sobj

#~4mins
endpoint_data = {}
for sname in sobj_names:
    print timestamp(sname)
    P = construct_poly(sname)
    A,B,C = get_slice_change_values(P)
    gapsA=[b-a for a,b in zip(A[1:-2], A[2:-1])]
    gapsB=[b-a for a,b in zip(B[1:-2], B[2:-1])]
    gapsC=[b-a for a,b in zip(C[1:-2], C[2:-1])]
    endpoint_data[sname] = (A,B,C,gapsA,gapsB,gapsC)

[2017-08-20 19:25:48.089313] sobjs/nitrosomonas.europaea.atcc.19718_BX321863_b.sobj [2017-08-20 19:25:50.326926] [get_slice_] Start [2017-08-20 19:25:50.415892] [get_slice_] Stop [2017-08-20 19:25:50.421448] sobjs/lactobacillus.casei.atcc.334_CP000423_j.sobj [2017-08-20 19:25:52.532934] [get_slice_] Start [2017-08-20 19:25:52.619583] [get_slice_] Stop [2017-08-20 19:25:52.622666] sobjs/methanosarcina.barkeri.str..fusaro_CP000099_bi.sobj [2017-08-20 19:26:07.299949] [get_slice_] Start [2017-08-20 19:26:07.463182] [get_slice_] Stop [2017-08-20 19:26:07.467222] sobjs/uncultured.alpha.proteobacterium_AY933788.sobj [2017-08-20 19:26:10.331763] [get_slice_] Start [2017-08-20 19:26:10.467244] [get_slice_] Stop [2017-08-20 19:26:10.472395] sobjs/tropheryma.whipplei.str..twist_AE014184_h.sobj [2017-08-20 19:26:12.889820] [get_slice_] Start [2017-08-20 19:26:12.999317] [get_slice_] Stop [2017-08-20 19:26:13.002867] sobjs/uncultured.thiotrichales.bacterium_AY934254_a.sobj [2017-08-20 19:26:15.323149] [get_slice_] Start [2017-08-20 19:26:15.419185] [get_slice_] Stop [2017-08-20 19:26:15.422628] sobjs/synechococcus.sp..wh.8102_BX569691.sobj [2017-08-20 19:26:17.543067] [get_slice_] Start [2017-08-20 19:26:17.633452] [get_slice_] Stop [2017-08-20 19:26:17.637381] sobjs/yersinia.pestis.kim_AE009952_ae.sobj [2017-08-20 19:26:20.343124] [get_slice_] Start [2017-08-20 19:26:20.451253] [get_slice_] Stop [2017-08-20 19:26:20.455165] sobjs/uncultured.alpha.proteobacterium_DQ093144_a.sobj [2017-08-20 19:26:22.817330] [get_slice_] Start [2017-08-20 19:26:22.922164] [get_slice_] Stop [2017-08-20 19:26:22.925869] sobjs/carboxydothermus.hydrogenoformans.z-2901_CP000141_af.sobj [2017-08-20 19:27:22.758843] [get_slice_] Start [2017-08-20 19:27:22.955449] [get_slice_] Stop [2017-08-20 19:27:22.959559] sobjs/haloferax.volcanii_HALTRK1.sobj [2017-08-20 19:27:25.594010] [get_slice_] Start [2017-08-20 19:27:25.727046] [get_slice_] Stop [2017-08-20 19:27:25.730525] sobjs/homo.sapiens_X04779.sobj [2017-08-20 19:27:28.267133] [get_slice_] Start [2017-08-20 19:27:28.385827] [get_slice_] Stop [2017-08-20 19:27:28.389454] sobjs/arabidopsis.thaliana_AC006340.sobj [2017-08-20 19:27:30.672684] [get_slice_] Start [2017-08-20 19:27:30.776751] [get_slice_] Stop [2017-08-20 19:27:30.780374] sobjs/chlamydophila.pneumoniae.tw-183_AE017159_j.sobj [2017-08-20 19:27:32.761150] [get_slice_] Start [2017-08-20 19:27:32.854725] [get_slice_] Stop [2017-08-20 19:27:32.858768] sobjs/thermoplasma.volcanium.gss1_BA000011_b.sobj [2017-08-20 19:27:35.182134] [get_slice_] Start [2017-08-20 19:27:35.285598] [get_slice_] Stop [2017-08-20 19:27:35.289144] sobjs/uncultured.gamma.proteobacterium_AY934184_a.sobj [2017-08-20 19:28:05.312856] [get_slice_] Start [2017-08-20 19:28:05.470721] [get_slice_] Stop [2017-08-20 19:28:05.475794] sobjs/nocardia.farcinica.ifm.10152_AP006618_az.sobj [2017-08-20 19:28:08.462880] [get_slice_] Start [2017-08-20 19:28:08.623237] [get_slice_] Stop [2017-08-20 19:28:08.627139] sobjs/chlorobium.phaeobacteroides.dsm.266_CP000492_ab.sobj [2017-08-20 19:28:10.705390] [get_slice_] Start [2017-08-20 19:28:10.820647] [get_slice_] Stop [2017-08-20 19:28:10.824072] sobjs/scenedesmus.obliquus_DQ396875_n.sobj [2017-08-20 19:28:12.585643] [get_slice_] Start [2017-08-20 19:28:12.679634] [get_slice_] Stop [2017-08-20 19:28:12.682376] sobjs/pyrococcus.abyssi.ge5_AJ248288.sobj [2017-08-20 19:28:15.104235] [get_slice_] Start [2017-08-20 19:28:15.217110] [get_slice_] Stop [2017-08-20 19:28:15.220917] sobjs/euglena.longa_AJ294725_s.sobj [2017-08-20 19:28:17.215377] [get_slice_] Start [2017-08-20 19:28:17.311297] [get_slice_] Stop [2017-08-20 19:28:17.314411] sobjs/marchantia.polymorpha_M.CAU.C.Mar.polym.1.sobj [2017-08-20 19:28:19.864133] [get_slice_] Start [2017-08-20 19:28:19.963306] [get_slice_] Stop [2017-08-20 19:28:19.966599] sobjs/lactobacillus.delbrueckii.subsp..bulgaricus.atcc.baa-365_CP000412_cj.sobj [2017-08-20 19:28:34.490477] [get_slice_] Start [2017-08-20 19:28:34.660189] [get_slice_] Stop [2017-08-20 19:28:34.663920] sobjs/uncultured.alteromonadales.bacterium_AY934393.sobj [2017-08-20 19:28:37.521662] [get_slice_] Start [2017-08-20 19:28:37.671556] [get_slice_] Stop [2017-08-20 19:28:37.675288] sobjs/methanospirillum.hungatei.jf-1_CP000254_ad.sobj [2017-08-20 19:28:39.999724] [get_slice_] Start [2017-08-20 19:28:40.124058] [get_slice_] Stop [2017-08-20 19:28:40.128975] sobjs/drosophila.melanogaster_DROTRF2.sobj [2017-08-20 19:28:42.839023] [get_slice_] Start [2017-08-20 19:28:42.978451] [get_slice_] Stop [2017-08-20 19:28:42.982621] sobjs/pelobacter.carbinolicus.dsm.2380_CP000142_af.sobj [2017-08-20 19:28:45.586979] [get_slice_] Start [2017-08-20 19:28:45.704073] [get_slice_] Stop [2017-08-20 19:28:45.708354] sobjs/bordetella.parapertussis_BX640433_e.sobj [2017-08-20 19:28:48.092333] [get_slice_] Start [2017-08-20 19:28:48.210402] [get_slice_] Stop [2017-08-20 19:28:48.214341] sobjs/magnetococcus.sp..mc-1_CP000471_ao.sobj [2017-08-20 19:28:50.597815] [get_slice_] Start [2017-08-20 19:28:50.704765] [get_slice_] Stop [2017-08-20 19:28:50.708363] sobjs/yarrowia.lipolytica.clib122_CR382129_by.sobj [2017-08-20 19:29:06.233334] [get_slice_] Start [2017-08-20 19:29:06.372483] [get_slice_] Stop [2017-08-20 19:29:06.375894] sobjs/hyperthermus.butylicus.dsm.5456_CP000493_ad.sobj [2017-08-20 19:29:09.087825] [get_slice_] Start [2017-08-20 19:29:09.232139] [get_slice_] Stop [2017-08-20 19:29:09.237952] sobjs/homo.sapiens_AC004932_k.sobj [2017-08-20 19:29:12.084545] [get_slice_] Start [2017-08-20 19:29:12.216894] [get_slice_] Stop [2017-08-20 19:29:12.220819] sobjs/uncultured.gamma.proteobacterium_AY934387_a.sobj [2017-08-20 19:29:14.727799] [get_slice_] Start [2017-08-20 19:29:14.848511] [get_slice_] Stop [2017-08-20 19:29:14.852431] sobjs/rickettsia.conorii.str..malish.7_AE008623.sobj [2017-08-20 19:29:17.048417] [get_slice_] Start [2017-08-20 19:29:17.156456] [get_slice_] Stop [2017-08-20 19:29:17.160064] sobjs/solibacter.usitatus.ellin6076_CP000473_ay.sobj [2017-08-20 19:29:19.508200] [get_slice_] Start [2017-08-20 19:29:19.622634] [get_slice_] Stop [2017-08-20 19:29:19.627277] sobjs/wigglesworthia.glossinidia.endosymbiont.of.glossina.brevipalpis_BA000021_ad.sobj [2017-08-20 19:29:34.570975] [get_slice_] Start [2017-08-20 19:29:34.733406] [get_slice_] Stop [2017-08-20 19:29:34.737772] sobjs/sulfolobus.tokodaii.str..7_BA000023_ak.sobj [2017-08-20 19:29:37.203934] [get_slice_] Start [2017-08-20 19:29:37.341358] [get_slice_] Stop [2017-08-20 19:29:37.346363] sobjs/uncultured.marinobacter.sp._AY934351.sobj [2017-08-20 19:29:40.153155] [get_slice_] Start [2017-08-20 19:29:40.302007] [get_slice_] Stop [2017-08-20 19:29:40.306011] sobjs/nephroselmis.olivacea_AF137379_g.sobj [2017-08-20 19:29:42.309983] [get_slice_] Start [2017-08-20 19:29:42.425109] [get_slice_] Stop

nacuts=[]
nbcuts=[]
nccuts=[]
for n in endpoint_data:
    a,b,c,ga,gb,gc = endpoint_data[n]
    nacuts.append(len(a))
    nbcuts.append(len(b))
    nccuts.append(len(c))
print "avg nacuts = {0}".format(avg(nacuts))
print "avg nbcuts = {0}".format(avg(nbcuts))
print "avg nccuts = {0}".format(avg(nccuts))

avg nacuts = 529.615384615 avg nbcuts = 200.564102564 avg nccuts = 367.897435897

nacuts

[568, 502, 526, 538, 545, 568, 588, 474, 549, 567, 472, 534, 546, 520, 597, 471, 439, 529, 499, 601, 532, 648, 484, 503, 555, 434, 505, 535, 524, 535, 547, 525, 488, 489, 497, 526, 586, 603, 506]

nd1={}
for sname in sobj_names:
    P = construct_poly(sname)
    nd1[sname] = len(P.d1_slices())

for n in nd1:
    a,b,c,ga,gb,gc = endpoint_data[n]
    print "{0}, {1}, {2}, {3}: [{4}]".format(len(a),len(b),len(c),nd1[n],n)

568, 209, 388, 534: [sobjs/methanosarcina.barkeri.str..fusaro_CP000099_bi.sobj] 502, 158, 318, 455: [sobjs/chlorobium.phaeobacteroides.dsm.266_CP000492_ab.sobj] 526, 185, 358, 501: [sobjs/homo.sapiens_X04779.sobj] 538, 221, 379, 535: [sobjs/solibacter.usitatus.ellin6076_CP000473_ay.sobj] 545, 217, 368, 566: [sobjs/uncultured.gamma.proteobacterium_AY934387_a.sobj] 568, 179, 381, 570: [sobjs/pyrococcus.abyssi.ge5_AJ248288.sobj] 588, 177, 402, 543: [sobjs/yersinia.pestis.kim_AE009952_ae.sobj] 474, 197, 347, 439: [sobjs/yarrowia.lipolytica.clib122_CR382129_by.sobj] 549, 228, 392, 519: [sobjs/sulfolobus.tokodaii.str..7_BA000023_ak.sobj] 567, 214, 393, 552: [sobjs/homo.sapiens_AC004932_k.sobj] 472, 177, 343, 487: [sobjs/euglena.longa_AJ294725_s.sobj] 534, 202, 374, 529: [sobjs/bordetella.parapertussis_BX640433_e.sobj] 546, 210, 382, 513: [sobjs/nitrosomonas.europaea.atcc.19718_BX321863_b.sobj] 520, 187, 370, 470: [sobjs/rickettsia.conorii.str..malish.7_AE008623.sobj] 597, 226, 426, 603: [sobjs/uncultured.alpha.proteobacterium_AY933788.sobj] 471, 212, 350, 453: [sobjs/methanospirillum.hungatei.jf-1_CP000254_ad.sobj] 439, 202, 328, 453: [sobjs/wigglesworthia.glossinidia.endosymbiont.of.glossina.brevipalpis_BA000021_ad.sobj] 529, 204, 361, 503: [sobjs/arabidopsis.thaliana_AC006340.sobj] 499, 159, 360, 507: [sobjs/synechococcus.sp..wh.8102_BX569691.sobj] 601, 201, 417, 553: [sobjs/drosophila.melanogaster_DROTRF2.sobj] 532, 226, 347, 558: [sobjs/magnetococcus.sp..mc-1_CP000471_ao.sobj] 648, 249, 420, 578: [sobjs/pelobacter.carbinolicus.dsm.2380_CP000142_af.sobj] 484, 187, 333, 469: [sobjs/nephroselmis.olivacea_AF137379_g.sobj] 503, 210, 368, 519: [sobjs/tropheryma.whipplei.str..twist_AE014184_h.sobj] 555, 213, 377, 541: [sobjs/carboxydothermus.hydrogenoformans.z-2901_CP000141_af.sobj] 434, 174, 305, 420: [sobjs/scenedesmus.obliquus_DQ396875_n.sobj] 505, 198, 379, 502: [sobjs/uncultured.gamma.proteobacterium_AY934184_a.sobj] 535, 192, 342, 494: [sobjs/lactobacillus.casei.atcc.334_CP000423_j.sobj] 524, 183, 346, 542: [sobjs/uncultured.alteromonadales.bacterium_AY934393.sobj] 535, 217, 418, 479: [sobjs/uncultured.thiotrichales.bacterium_AY934254_a.sobj] 547, 223, 400, 563: [sobjs/uncultured.marinobacter.sp._AY934351.sobj] 525, 211, 348, 504: [sobjs/uncultured.alpha.proteobacterium_DQ093144_a.sobj] 488, 184, 338, 449: [sobjs/chlamydophila.pneumoniae.tw-183_AE017159_j.sobj] 489, 198, 336, 502: [sobjs/haloferax.volcanii_HALTRK1.sobj] 497, 188, 342, 521: [sobjs/thermoplasma.volcanium.gss1_BA000011_b.sobj] 526, 181, 379, 534: [sobjs/marchantia.polymorpha_M.CAU.C.Mar.polym.1.sobj] 586, 194, 374, 600: [sobjs/nocardia.farcinica.ifm.10152_AP006618_az.sobj] 603, 229, 425, 582: [sobjs/hyperthermus.butylicus.dsm.5456_CP000493_ad.sobj] 506, 200, 334, 482: [sobjs/lactobacillus.delbrueckii.subsp..bulgaricus.atcc.baa-365_CP000412_cj.sobj]

minas=[]
maxas=[]
minbs=[]
maxbs=[]
mincs=[]
maxcs=[]
for n in endpoint_data:
    a,b,c,ga,gb,gc = endpoint_data[n]
    maxas.append(a[-2])
    minas.append(a[1])
    maxbs.append(b[-2])
    minbs.append(b[1])
    maxcs.append(c[-2])
    mincs.append(c[1])

print "Avg max as = {0}".format(avg(maxas))
print "Std max as = {0}".format(std(maxas))
print "Avg min as = {0}".format(avg(minas))
print "Std min as = {0}".format(std(minas))
print '*'*40
print "Avg max bs = {0}".format(avg(maxbs))
print "Std max bs = {0}".format(std(maxbs))
print "Avg min bs = {0}".format(avg(minbs))
print "Std min bs = {0}".format(std(minbs))
print '*'*40
print "Avg max cs = {0}".format(avg(maxcs))
print "Std max cs = {0}".format(std(maxcs))
print "Avg min cs = {0}".format(avg(mincs))
print "Std min cs = {0}".format(std(mincs))

Avg max as = 346.742564103 Std max as = 83.6529514773 Avg min as = -5.38763751843 Std min as = 6.4524371905 **************************************** Avg max bs = 5.46581196581 Std max bs = 2.09279713286 Avg min bs = -4.67435897436 Std min bs = 1.58213875833 **************************************** Avg max cs = 1.80460357205 Std max cs = 2.37541030917 Avg min cs = -38.6484615385 Std min cs = 9.83316902213

posl =  [x-1 for x in positive_orthant_slices.values()]
print posl
print avg(posl)
print std(posl)

[2, 4, 5, 4, 5, 8, 3, 9, 3, 4, 4, 8, 5, 3, 6, 3, 5, 8, 7, 5, 5, 4, 6, 2, 7, 11, 3, 6, 3, 3, 3, 4, 3, 2, 2, 4, 3, 7, 9] 4.82051282051 2.21731871888

%md
# On the number of changes of slices.
```Avg max as = 346.742564103
Std max as = 83.6529514773
Avg min as = -5.38763751843
Std min as = 6.4524371905
****************************************
Avg max bs = 5.46581196581
Std max bs = 2.09279713286
Avg min bs = -4.67435897436
Std min bs = 1.58213875833
****************************************
Avg max cs = 1.80460357205
Std max cs = 2.37541030917
Avg min cs = -38.6484615385
Std min cs = 9.83316902213```
# On the the number of slices that intersect the positive orthant
```[2, 4, 5, 4, 5, 8, 3, 9, 3, 4, 4, 8, 5, 3, 6, 3, 5, 8, 7, 5, 5, 4, 6, 2, 7, 11, 3, 6, 3, 3, 3, 4, 3, 2, 2, 4, 3, 7, 9]
4.82051282051
2.21731871888```

On the number of changes of slices.

Std max as = 83.6529514773
Avg min as = -5.38763751843
Std min as = 6.4524371905
****************************************
Avg max bs = 5.46581196581
Std max bs = 2.09279713286
Avg min bs = -4.67435897436
Std min bs = 1.58213875833
****************************************
Avg max cs = 1.80460357205
Std max cs = 2.37541030917
Avg min cs = -38.6484615385
Std min cs = 9.83316902213```
# On the the number of slices that intersect the positive orthant
```[2, 4, 5, 4, 5, 8, 3, 9, 3, 4, 4, 8, 5, 3, 6, 3, 5, 8, 7, 5, 5, 4, 6, 2, 7, 11, 3, 6, 3, 3, 3, 4, 3, 2, 2, 4, 3, 7, 9]
4.82051282051
2.21731871888```

min(minas)

-118/5

avg(A[1:-1])
std(A[1:-1])
avg(gapsA)
std(gapsA)

41.715657207017451 47.63756932504068 0.68121546961325974 3.4700515830513914

avg(B[1:-1])
std(B[1:-1])
avg(gapsB)
std(gapsB)

-0.69641424037563737 1.8270376336988228 0.049323671497584543 0.082031006522740951

avg(C[1:-1])
std(C[1:-1])
avg(gapsC)
std(gapsC)

-14.457445259293465 8.1577212435387221 0.12262532981530345 0.34094968993165325

A,B,C = get_slice_change_values(P)

avg(A[1:-1])
std(A[1:-1])

41.715657207017451 47.63756932504068

avg(B[1:-1])
std(B[1:-1])

-0.41780050745608016 1.5332737683607829

avg(C[1:-1])
std(C[1:-1])

-13.548168507554537 6.164310192802259

criteria = lambda sl: is_bounded(sl) and slice_has_positive_slope(sl)
slb=None
for snameb in sobj_names:
    P=construct_poly(snameb)
    slb = find_slice(P,criteria)
    if slb:
        break

#print slb.original_vertex
#print snameb
#slb.show()

#plot_slice_cuts(slb)

positive_orthant = Polyhedron(ieqs=[[0, 1, 0, 0], [0, 0, 1, 0], [0, 0, 0, 1]],base_ring=QQ)

#positive_orthant.show()

print sl.original_vertex

A vertex at (2, 0, 8, -147/10)

sname

'sobjs/lactobacillus.casei.atcc.334_CP000423_j.sobj'

#sl.show()

#plot_slice_cuts(sl)

type(sl.rays_list()[0])

<type 'list'>

print sobj_names[:20]
print '*'*40
print sobj_names[20:]

['sobjs/nitrosomonas.europaea.atcc.19718_BX321863_b.sobj', 'sobjs/lactobacillus.casei.atcc.334_CP000423_j.sobj', 'sobjs/methanosarcina.barkeri.str..fusaro_CP000099_bi.sobj', 'sobjs/uncultured.alpha.proteobacterium_AY933788.sobj', 'sobjs/tropheryma.whipplei.str..twist_AE014184_h.sobj', 'sobjs/uncultured.thiotrichales.bacterium_AY934254_a.sobj', 'sobjs/synechococcus.sp..wh.8102_BX569691.sobj', 'sobjs/yersinia.pestis.kim_AE009952_ae.sobj', 'sobjs/uncultured.alpha.proteobacterium_DQ093144_a.sobj', 'sobjs/carboxydothermus.hydrogenoformans.z-2901_CP000141_af.sobj', 'sobjs/haloferax.volcanii_HALTRK1.sobj', 'sobjs/homo.sapiens_X04779.sobj', 'sobjs/arabidopsis.thaliana_AC006340.sobj', 'sobjs/chlamydophila.pneumoniae.tw-183_AE017159_j.sobj', 'sobjs/thermoplasma.volcanium.gss1_BA000011_b.sobj', 'sobjs/uncultured.gamma.proteobacterium_AY934184_a.sobj', 'sobjs/nocardia.farcinica.ifm.10152_AP006618_az.sobj', 'sobjs/chlorobium.phaeobacteroides.dsm.266_CP000492_ab.sobj', 'sobjs/scenedesmus.obliquus_DQ396875_n.sobj', 'sobjs/pyrococcus.abyssi.ge5_AJ248288.sobj'] **************************************** ['sobjs/euglena.longa_AJ294725_s.sobj', 'sobjs/marchantia.polymorpha_M.CAU.C.Mar.polym.1.sobj', 'sobjs/lactobacillus.delbrueckii.subsp..bulgaricus.atcc.baa-365_CP000412_cj.sobj', 'sobjs/uncultured.alteromonadales.bacterium_AY934393.sobj', 'sobjs/methanospirillum.hungatei.jf-1_CP000254_ad.sobj', 'sobjs/drosophila.melanogaster_DROTRF2.sobj', 'sobjs/pelobacter.carbinolicus.dsm.2380_CP000142_af.sobj', 'sobjs/bordetella.parapertussis_BX640433_e.sobj', 'sobjs/magnetococcus.sp..mc-1_CP000471_ao.sobj', 'sobjs/yarrowia.lipolytica.clib122_CR382129_by.sobj', 'sobjs/hyperthermus.butylicus.dsm.5456_CP000493_ad.sobj', 'sobjs/homo.sapiens_AC004932_k.sobj', 'sobjs/uncultured.gamma.proteobacterium_AY934387_a.sobj', 'sobjs/rickettsia.conorii.str..malish.7_AE008623.sobj', 'sobjs/solibacter.usitatus.ellin6076_CP000473_ay.sobj', 'sobjs/wigglesworthia.glossinidia.endosymbiont.of.glossina.brevipalpis_BA000021_ad.sobj', 'sobjs/sulfolobus.tokodaii.str..7_BA000023_ak.sobj', 'sobjs/uncultured.marinobacter.sp._AY934351.sobj', 'sobjs/nephroselmis.olivacea_AF137379_g.sobj']

# ~10 min
pos_slope_regions = {}
for sname in sobj_names[:20]:
    P=construct_poly(sname)
    print timestamp(sname)
    pos_slope_regions[sname] = collect_positive_slope_regions(P)
print "Done!"

[2017-08-19 23:38:41.606441] sobjs/nitrosomonas.europaea.atcc.19718_BX321863_b.sobj [2017-08-19 23:38:41.614410] [collect_po] Start [2017-08-19 23:39:03.716789] [collect_po] Stop [2017-08-19 23:39:07.152359] sobjs/lactobacillus.casei.atcc.334_CP000423_j.sobj [2017-08-19 23:39:07.160310] [collect_po] Start [2017-08-19 23:39:26.086120] [collect_po] Stop [2017-08-19 23:39:30.366061] sobjs/methanosarcina.barkeri.str..fusaro_CP000099_bi.sobj [2017-08-19 23:39:30.375563] [collect_po] Start [2017-08-19 23:39:53.557541] [collect_po] Stop [2017-08-19 23:39:57.153297] sobjs/uncultured.alpha.proteobacterium_AY933788.sobj [2017-08-19 23:39:57.165216] [collect_po] Start [2017-08-19 23:40:27.773250] [collect_po] Stop [2017-08-19 23:40:31.232598] sobjs/tropheryma.whipplei.str..twist_AE014184_h.sobj [2017-08-19 23:40:31.240080] [collect_po] Start [2017-08-19 23:40:55.064928] [collect_po] Stop [2017-08-19 23:40:58.741846] sobjs/uncultured.thiotrichales.bacterium_AY934254_a.sobj [2017-08-19 23:40:58.750378] [collect_po] Start [2017-08-19 23:41:19.174682] [collect_po] Stop [2017-08-19 23:41:23.196183] sobjs/synechococcus.sp..wh.8102_BX569691.sobj [2017-08-19 23:41:23.202453] [collect_po] Start [2017-08-19 23:41:40.347920] [collect_po] Stop [2017-08-19 23:41:45.273100] sobjs/yersinia.pestis.kim_AE009952_ae.sobj [2017-08-19 23:41:45.285925] [collect_po] Start A vertex at (8, 0, 25, 6151/100) ((((((...)((....)))(((....)(...))(((...)(...))(...)))))(...))(...)(...)) [2017-08-19 23:42:10.086407] [collect_po] Stop [2017-08-19 23:42:14.751761] sobjs/uncultured.alpha.proteobacterium_DQ093144_a.sobj [2017-08-19 23:42:14.762471] [collect_po] Start A vertex at (5, 4, 20, 6711/100) ((...)((><(...)((...)><(...)((...)(....)))(...)(...))(...)(...)(.....)))> [2017-08-19 23:42:39.697426] [collect_po] Stop [2017-08-19 23:42:44.681820] sobjs/carboxydothermus.hydrogenoformans.z-2901_CP000141_af.sobj [2017-08-19 23:42:44.698942] [collect_po] Start [2017-08-19 23:43:12.148576] [collect_po] Stop [2017-08-19 23:43:17.369569] sobjs/haloferax.volcanii_HALTRK1.sobj [2017-08-19 23:43:17.382424] [collect_po] Start [2017-08-19 23:43:43.041361] [collect_po] Stop [2017-08-19 23:43:44.926001] sobjs/homo.sapiens_X04779.sobj [2017-08-19 23:43:44.931955] [collect_po] Start A vertex at (6, 5, 22, 66) ((...)>(>((...)(((>(...)(...)>(...)>(...)(...))(...))(...)(...)))(...)))>. [2017-08-19 23:44:17.443214] [collect_po] Stop [2017-08-19 23:44:19.371956] sobjs/arabidopsis.thaliana_AC006340.sobj [2017-08-19 23:44:19.379875] [collect_po] Start A vertex at (6, 0, 22, 61) ((...)(...)(...)(...)((((....)(((...)(....))(...))((...)(....)))(...))))> [2017-08-19 23:44:51.764457] [collect_po] Stop [2017-08-19 23:44:53.377878] sobjs/chlamydophila.pneumoniae.tw-183_AE017159_j.sobj [2017-08-19 23:44:53.385943] [collect_po] Start [2017-08-19 23:45:19.953656] [collect_po] Stop [2017-08-19 23:45:25.823992] sobjs/thermoplasma.volcanium.gss1_BA000011_b.sobj [2017-08-19 23:45:25.834406] [collect_po] Start A vertex at (8, 1, 25, 5721/100) (((((...)((...)((...)(((...)(...)(.....))(....)))))(...))>(...))((...)))> A vertex at (6, 5, 22, 623/10) ((...)(...)>(>(((...)><(...)(...))(...)((...)(...)((...)(...)>)))(...)))> [2017-08-19 23:45:57.790974] [collect_po] Stop [2017-08-19 23:45:59.553567] sobjs/uncultured.gamma.proteobacterium_AY934184_a.sobj [2017-08-19 23:45:59.560183] [collect_po] Start [2017-08-19 23:46:33.348490] [collect_po] Stop [2017-08-19 23:46:35.698854] sobjs/nocardia.farcinica.ifm.10152_AP006618_az.sobj [2017-08-19 23:46:35.708975] [collect_po] Start [2017-08-19 23:47:11.259093] [collect_po] Stop [2017-08-19 23:47:17.709421] sobjs/chlorobium.phaeobacteroides.dsm.266_CP000492_ab.sobj [2017-08-19 23:47:17.715382] [collect_po] Start [2017-08-19 23:47:42.087884] [collect_po] Stop [2017-08-19 23:47:43.420894] sobjs/scenedesmus.obliquus_DQ396875_n.sobj [2017-08-19 23:47:43.426226] [collect_po] Start [2017-08-19 23:48:03.920495] [collect_po] Stop [2017-08-19 23:48:06.068341] sobjs/pyrococcus.abyssi.ge5_AJ248288.sobj [2017-08-19 23:48:06.074504] [collect_po] Start [2017-08-19 23:48:45.438665] [collect_po] Stop Done!

# ~16 min
for sname in sobj_names[20:]:
    P=construct_poly(sname)
    print timestamp(sname)
    pos_slope_regions[sname] = collect_positive_slope_regions(P)
print "Done!"

[2017-08-19 23:56:51.176019] sobjs/euglena.longa_AJ294725_s.sobj [2017-08-19 23:56:51.183036] [collect_po] Start [2017-08-19 23:57:21.307958] [collect_po] Stop [2017-08-19 23:57:23.171075] sobjs/marchantia.polymorpha_M.CAU.C.Mar.polym.1.sobj [2017-08-19 23:57:23.176998] [collect_po] Start [2017-08-19 23:58:01.488527] [collect_po] Stop [2017-08-19 23:58:03.681874] sobjs/lactobacillus.delbrueckii.subsp..bulgaricus.atcc.baa-365_CP000412_cj.sobj [2017-08-19 23:58:03.692124] [collect_po] Start [2017-08-19 23:58:34.335928] [collect_po] Stop [2017-08-19 23:58:36.173759] sobjs/uncultured.alteromonadales.bacterium_AY934393.sobj [2017-08-19 23:58:36.180933] [collect_po] Start [2017-08-19 23:59:17.077055] [collect_po] Stop [2017-08-19 23:59:18.784276] sobjs/methanospirillum.hungatei.jf-1_CP000254_ad.sobj [2017-08-19 23:59:18.794376] [collect_po] Start [2017-08-19 23:59:52.590177] [collect_po] Stop [2017-08-19 23:59:54.681953] sobjs/drosophila.melanogaster_DROTRF2.sobj [2017-08-19 23:59:54.688484] [collect_po] Start [2017-08-20 00:00:33.241219] [collect_po] Stop [2017-08-20 00:00:42.518491] sobjs/pelobacter.carbinolicus.dsm.2380_CP000142_af.sobj [2017-08-20 00:00:42.529357] [collect_po] Start [2017-08-20 00:01:27.016357] [collect_po] Stop [2017-08-20 00:01:28.889021] sobjs/bordetella.parapertussis_BX640433_e.sobj [2017-08-20 00:01:28.895889] [collect_po] Start [2017-08-20 00:02:11.992077] [collect_po] Stop [2017-08-20 00:02:21.368866] sobjs/magnetococcus.sp..mc-1_CP000471_ao.sobj [2017-08-20 00:02:21.380157] [collect_po] Start [2017-08-20 00:03:06.134501] [collect_po] Stop [2017-08-20 00:03:07.789930] sobjs/yarrowia.lipolytica.clib122_CR382129_by.sobj [2017-08-20 00:03:07.797203] [collect_po] Start A vertex at (7, 0, 23, 621/10) (((...)((((.....)((...)(....)))(...))(...))((...)(...)))(.......)(....))> [2017-08-20 00:03:42.871498] [collect_po] Stop [2017-08-20 00:03:44.754622] sobjs/hyperthermus.butylicus.dsm.5456_CP000493_ad.sobj [2017-08-20 00:03:44.762735] [collect_po] Start [2017-08-20 00:04:31.767547] [collect_po] Stop [2017-08-20 00:04:42.034860] sobjs/homo.sapiens_AC004932_k.sobj [2017-08-20 00:04:42.041925] [collect_po] Start [2017-08-20 00:05:25.992114] [collect_po] Stop [2017-08-20 00:05:28.246808] sobjs/uncultured.gamma.proteobacterium_AY934387_a.sobj [2017-08-20 00:05:28.255871] [collect_po] Start [2017-08-20 00:06:16.377639] [collect_po] Stop [2017-08-20 00:06:18.117221] sobjs/rickettsia.conorii.str..malish.7_AE008623.sobj [2017-08-20 00:06:18.124040] [collect_po] Start [2017-08-20 00:07:12.544709] [collect_po] Stop [2017-08-20 00:07:15.273720] sobjs/solibacter.usitatus.ellin6076_CP000473_ay.sobj [2017-08-20 00:07:15.287647] [collect_po] Start [2017-08-20 00:08:03.219470] [collect_po] Stop [2017-08-20 00:08:04.790433] sobjs/wigglesworthia.glossinidia.endosymbiont.of.glossina.brevipalpis_BA000021_ad.sobj [2017-08-20 00:08:04.800037] [collect_po] Start [2017-08-20 00:08:44.180285] [collect_po] Stop [2017-08-20 00:08:45.909214] sobjs/sulfolobus.tokodaii.str..7_BA000023_ak.sobj [2017-08-20 00:08:45.918954] [collect_po] Start [2017-08-20 00:09:31.996786] [collect_po] Stop [2017-08-20 00:09:43.768470] sobjs/uncultured.marinobacter.sp._AY934351.sobj [2017-08-20 00:09:43.778739] [collect_po] Start [2017-08-20 00:10:40.436942] [collect_po] Stop [2017-08-20 00:10:41.943840] sobjs/nephroselmis.olivacea_AF137379_g.sobj [2017-08-20 00:10:41.951010] [collect_po] Start A vertex at (7, 0, 23, 294/5) .<(((((...)(...))((....)(....)))((...)(....)(...)((...)(....))))(....))> [2017-08-20 00:11:20.587018] [collect_po] Stop Done!

positive_slope_pcts = []
for n in pos_slope_regions:
    nr,nb,nbp,nu,nup = pos_slope_regions[n]
    #print "{0}: {1}; {2}/{3}".format(n,pos_slope_regions[n],nup+nbp,nr)
    print "{0}/{1}={2}%".format(nup+nbp,nr,100.0*float(nup+nbp)/nr)
    positive_slope_pcts.append(100.0*float(nup+nbp)/nr)

21/534=3.93258426966292 9/455=1.97802197802198 15/501=2.99401197604790 4/535=0.747663551401869 4/566=0.706713780918728 16/570=2.80701754385965 18/543=3.31491712707182 15/439=3.41685649202733 16/519=3.08285163776493 19/552=3.44202898550725 13/487=2.66940451745380 8/529=1.51228733459357 0/513=0.000000000000000 7/470=1.48936170212766 5/603=0.829187396351575 2/453=0.441501103752759 15/453=3.31125827814570 6/503=1.19284294234592 11/507=2.16962524654832 22/553=3.97830018083183 6/558=1.07526881720430 13/578=2.24913494809689 8/469=1.70575692963753 14/519=2.69749518304432 9/541=1.66358595194085 7/420=1.66666666666667 0/502=0.000000000000000 7/494=1.41700404858300 6/542=1.10701107011070 27/479=5.63674321503132 6/563=1.06571936056838 12/504=2.38095238095238 13/449=2.89532293986637 20/502=3.98406374501992 4/521=0.767754318618042 27/534=5.05617977528090 5/600=0.833333333333333 11/582=1.89003436426117 8/482=1.65975103734440

positive_slope_pcts
avg(positive_slope_pcts)
std(positive_slope_pcts)

[3.93258426966292, 1.97802197802198, 2.99401197604790, 0.747663551401869, 0.706713780918728, 2.80701754385965, 3.31491712707182, 3.41685649202733, 3.08285163776493, 3.44202898550725, 2.66940451745380, 1.51228733459357, 0.000000000000000, 1.48936170212766, 0.829187396351575, 0.441501103752759, 3.31125827814570, 1.19284294234592, 2.16962524654832, 3.97830018083183, 1.07526881720430, 2.24913494809689, 1.70575692963753, 2.69749518304432, 1.66358595194085, 1.66666666666667, 0.000000000000000, 1.41700404858300, 1.10701107011070, 5.63674321503132, 1.06571936056838, 2.38095238095238, 2.89532293986637, 3.98406374501992, 0.767754318618042, 5.05617977528090, 0.833333333333333, 1.89003436426117, 1.65975103734440] 2.1479029264101532 1.3203760095624049

%md
# Percentage of positive slope slices
```[3.93258426966292, 1.97802197802198, 2.99401197604790, 0.747663551401869, 0.706713780918728, 2.80701754385965, 3.31491712707182, 3.41685649202733, 3.08285163776493, 3.44202898550725, 2.66940451745380, 1.51228733459357, 0.000000000000000, 1.48936170212766, 0.829187396351575, 0.441501103752759, 3.31125827814570, 1.19284294234592, 2.16962524654832, 3.97830018083183, 1.07526881720430, 2.24913494809689, 1.70575692963753, 2.69749518304432, 1.66358595194085, 1.66666666666667, 0.000000000000000, 1.41700404858300, 1.10701107011070, 5.63674321503132, 1.06571936056838, 2.38095238095238, 2.89532293986637, 3.98406374501992, 0.767754318618042, 5.05617977528090, 0.833333333333333, 1.89003436426117, 1.65975103734440]
2.1479029264101532
1.3203760095624049```

Percentage of positive slope slices

2.1479029264101532
1.3203760095624049```

for n in pos_slope_regions:
    nr,nb,nbp,bp,up = pos_slope_regions[n]
    if len(up)>0:
        print n

for r in sl.ray_generator():
    break

r.vector().norm(1)

17

from itertools import combinations
distance = lambda u,v: sum([abs(x-y) for x,y in zip(u,v)])
diameter = lambda sl: max([distance(u,v) for u,v in combinations(sl.vertices_list(),2)])
diameter(sl)

57/10

binomial(16,2)

120

b_vals = [-10,-5,-2,-1,0,1,2,5,10]
acsp={}
print timestamp()
for sname in sobj_names:
    P = RNAPolytope.construct_from_file(sname)
    acsp[sname] = {'a':[],'w':[],'s':[],'p':[]}
    for b in b_vals:
        a,w,s,p=get_cones_stripes_polygons(P,b,signatures_only=True)
        acsp[sname]['a'].extend(a)
        acsp[sname]['w'].extend(w)
        acsp[sname]['s'].extend(s)
        acsp[sname]['p'].extend(p)
    acsp[sname]['ta'] = len(acsp[sname]['a'])
    acsp[sname]['tw'] = len(acsp[sname]['w'])
    acsp[sname]['ts'] = len(acsp[sname]['s'])
    acsp[sname]['tp'] = len(acsp[sname]['p'])
    acsp[sname]['tu'] = acsp[sname]['tw']+acsp[sname]['ts']
print timestamp()
# 12 minutes

[2017-08-10 21:50:42.416642] [2017-08-10 21:58:12.200428]

acsp_pcts ={}
keys = ['tp','tu','tw','ts']
for sname in acsp:
    acsp_pcts[sname] = { k[-1]:100.0*acsp[sname][k]/float(acsp[sname]['ta']) for k in keys}

acsp_pcts

{'sobjs/chlorobium.phaeobacteroides.dsm.266_CP000492_ab.sobj': {'p': 66.1676646706587, 's': 20.3592814371257, 'u': 33.8323353293413, 'w': 13.4730538922156}, 'sobjs/homo.sapiens_X04779.sobj': {'p': 67.2823218997361, 's': 20.8443271767810, 'u': 32.7176781002639, 'w': 11.8733509234828}, 'sobjs/solibacter.usitatus.ellin6076_CP000473_ay.sobj': {'p': 60.0000000000000, 's': 30.1369863013699, 'u': 40.0000000000000, 'w': 9.86301369863014}, 'sobjs/uncultured.gamma.proteobacterium_AY934387_a.sobj': {'p': 58.8235294117647, 's': 31.5508021390374, 'u': 41.1764705882353, 'w': 9.62566844919786}, 'sobjs/pyrococcus.abyssi.ge5_AJ248288.sobj': {'p': 70.1204819277108, 's': 19.0361445783133, 'u': 29.8795180722892, 'w': 10.8433734939759}, 'sobjs/yersinia.pestis.kim_AE009952_ae.sobj': {'p': 68.3937823834197, 's': 17.6165803108808, 'u': 31.6062176165803, 'w': 13.9896373056995}, 'sobjs/yarrowia.lipolytica.clib122_CR382129_by.sobj': {'p': 67.8160919540230, 's': 16.6666666666667, 'u': 32.1839080459770, 'w': 15.5172413793103}, 'sobjs/sulfolobus.tokodaii.str..7_BA000023_ak.sobj': {'p': 62.3655913978495, 's': 25.5376344086022, 'u': 37.6344086021505, 'w': 12.0967741935484}, 'sobjs/euglena.longa_AJ294725_s.sobj': {'p': 62.2222222222222, 's': 25.2777777777778, 'u': 37.7777777777778, 'w': 12.5000000000000}, 'sobjs/bordetella.parapertussis_BX640433_e.sobj': {'p': 67.6691729323308, 's': 21.0526315789474, 'u': 32.3308270676692, 'w': 11.2781954887218}, 'sobjs/nitrosomonas.europaea.atcc.19718_BX321863_b.sobj': {'p': 62.5352112676056, 's': 27.3239436619718, 'u': 37.4647887323944, 'w': 10.1408450704225}, 'sobjs/rickettsia.conorii.str..malish.7_AE008623.sobj': {'p': 68.8442211055276, 's': 22.1105527638191, 'u': 31.1557788944724, 'w': 9.04522613065327}, 'sobjs/uncultured.alpha.proteobacterium_AY933788.sobj': {'p': 66.0975609756098, 's': 22.9268292682927, 'u': 33.9024390243902, 'w': 10.9756097560976}, 'sobjs/methanospirillum.hungatei.jf-1_CP000254_ad.sobj': {'p': 62.5000000000000, 's': 24.1071428571429, 'u': 37.5000000000000, 'w': 13.3928571428571}, 'sobjs/wigglesworthia.glossinidia.endosymbiont.of.glossina.brevipalpis_BA000021_ad.sobj': {'p': 59.7014925373134, 's': 29.5522388059701, 'u': 40.2985074626866, 'w': 10.7462686567164}, 'sobjs/arabidopsis.thaliana_AC006340.sobj': {'p': 61.8384401114206, 's': 25.6267409470752, 'u': 38.1615598885794, 'w': 12.5348189415042}, 'sobjs/synechococcus.sp..wh.8102_BX569691.sobj': {'p': 56.7823343848580, 's': 31.8611987381703, 'u': 43.2176656151420, 'w': 11.3564668769716}, 'sobjs/uncultured.alteromonadales.bacterium_AY934393.sobj': {'p': 59.5441595441595, 's': 30.1994301994302, 'u': 40.4558404558405, 'w': 10.2564102564103}, 'sobjs/magnetococcus.sp..mc-1_CP000471_ao.sobj': {'p': 64.4155844155844, 's': 26.2337662337662, 'u': 35.5844155844156, 'w': 9.35064935064935}, 'sobjs/pelobacter.carbinolicus.dsm.2380_CP000142_af.sobj': {'p': 65.0366748166259, 's': 23.9608801955990, 'u': 34.9633251833741, 'w': 11.0024449877751}, 'sobjs/nephroselmis.olivacea_AF137379_g.sobj': {'p': 57.0532915360502, 's': 31.6614420062696, 'u': 42.9467084639498, 'w': 11.2852664576803}, 'sobjs/tropheryma.whipplei.str..twist_AE014184_h.sobj': {'p': 63.7055837563452, 's': 24.8730964467005, 'u': 36.2944162436548, 'w': 11.4213197969543}, 'sobjs/marchantia.polymorpha_M.CAU.C.Mar.polym.1.sobj': {'p': 68.0628272251309, 's': 17.8010471204188, 'u': 31.9371727748691, 'w': 14.1361256544503}, 'sobjs/homo.sapiens_AC004932_k.sobj': {'p': 63.4877384196185, 's': 24.2506811989101, 'u': 36.5122615803815, 'w': 12.2615803814714}, 'sobjs/scenedesmus.obliquus_DQ396875_n.sobj': {'p': 49.3975903614458, 's': 36.1445783132530, 'u': 50.6024096385542, 'w': 14.4578313253012}, 'sobjs/uncultured.gamma.proteobacterium_AY934184_a.sobj': {'p': 68.6111111111111, 's': 18.8888888888889, 'u': 31.3888888888889, 'w': 12.5000000000000}, 'sobjs/lactobacillus.casei.atcc.334_CP000423_j.sobj': {'p': 62.0498614958449, 's': 25.4847645429363, 'u': 37.9501385041551, 'w': 12.4653739612188}, 'sobjs/drosophila.melanogaster_DROTRF2.sobj': {'p': 65.6167979002625, 's': 24.9343832020997, 'u': 34.3832020997375, 'w': 9.44881889763779}, 'sobjs/carboxydothermus.hydrogenoformans.z-2901_CP000141_af.sobj': {'p': 63.3245382585752, 's': 27.1767810026385, 'u': 36.6754617414248, 'w': 9.49868073878628}, 'sobjs/uncultured.marinobacter.sp._AY934351.sobj': {'p': 66.7493796526055, 's': 24.3176178660050, 'u': 33.2506203473945, 'w': 8.93300248138958}, 'sobjs/uncultured.alpha.proteobacterium_DQ093144_a.sobj': {'p': 64.6067415730337, 's': 22.7528089887640, 'u': 35.3932584269663, 'w': 12.6404494382022}, 'sobjs/chlamydophila.pneumoniae.tw-183_AE017159_j.sobj': {'p': 69.6569920844327, 's': 18.4696569920844, 'u': 30.3430079155673, 'w': 11.8733509234828}, 'sobjs/haloferax.volcanii_HALTRK1.sobj': {'p': 63.6871508379888, 's': 23.7430167597765, 'u': 36.3128491620112, 'w': 12.5698324022346}, 'sobjs/thermoplasma.volcanium.gss1_BA000011_b.sobj': {'p': 61.7816091954023, 's': 25.2873563218391, 'u': 38.2183908045977, 'w': 12.9310344827586}, 'sobjs/uncultured.thiotrichales.bacterium_AY934254_a.sobj': {'p': 64.8044692737430, 's': 22.6256983240223, 'u': 35.1955307262570, 'w': 12.5698324022346}, 'sobjs/methanosarcina.barkeri.str..fusaro_CP000099_bi.sobj': {'p': 66.6666666666667, 's': 19.1601049868766, 'u': 33.3333333333333, 'w': 14.1732283464567}, 'sobjs/nocardia.farcinica.ifm.10152_AP006618_az.sobj': {'p': 66.9250645994832, 's': 21.4470284237726, 'u': 33.0749354005168, 'w': 11.6279069767442}, 'sobjs/hyperthermus.butylicus.dsm.5456_CP000493_ad.sobj': {'p': 62.2691292875989, 's': 28.2321899736148, 'u': 37.7308707124011, 'w': 9.49868073878628}, 'sobjs/lactobacillus.delbrueckii.subsp..bulgaricus.atcc.baa-365_CP000412_cj.sobj': {'p': 64.5604395604396, 's': 23.0769230769231, 'u': 35.4395604395604, 'w': 12.3626373626374}}

data = [[],[],[],[]]
for i,typ in enumerate(['p','u','s','w']):
    for sname in acsp_pcts:
        data[i].append(acsp_pcts[sname][typ])

data

[[66.1676646706587, 67.2823218997361, 60.0000000000000, 58.8235294117647, 70.1204819277108, 68.3937823834197, 67.8160919540230, 62.3655913978495, 62.2222222222222, 67.6691729323308, 62.5352112676056, 68.8442211055276, 66.0975609756098, 62.5000000000000, 59.7014925373134, 61.8384401114206, 56.7823343848580, 59.5441595441595, 64.4155844155844, 65.0366748166259, 57.0532915360502, 63.7055837563452, 68.0628272251309, 63.4877384196185, 49.3975903614458, 68.6111111111111, 62.0498614958449, 65.6167979002625, 63.3245382585752, 66.7493796526055, 64.6067415730337, 69.6569920844327, 63.6871508379888, 61.7816091954023, 64.8044692737430, 66.6666666666667, 66.9250645994832, 62.2691292875989, 64.5604395604396], [33.8323353293413, 32.7176781002639, 40.0000000000000, 41.1764705882353, 29.8795180722892, 31.6062176165803, 32.1839080459770, 37.6344086021505, 37.7777777777778, 32.3308270676692, 37.4647887323944, 31.1557788944724, 33.9024390243902, 37.5000000000000, 40.2985074626866, 38.1615598885794, 43.2176656151420, 40.4558404558405, 35.5844155844156, 34.9633251833741, 42.9467084639498, 36.2944162436548, 31.9371727748691, 36.5122615803815, 50.6024096385542, 31.3888888888889, 37.9501385041551, 34.3832020997375, 36.6754617414248, 33.2506203473945, 35.3932584269663, 30.3430079155673, 36.3128491620112, 38.2183908045977, 35.1955307262570, 33.3333333333333, 33.0749354005168, 37.7308707124011, 35.4395604395604], [20.3592814371257, 20.8443271767810, 30.1369863013699, 31.5508021390374, 19.0361445783133, 17.6165803108808, 16.6666666666667, 25.5376344086022, 25.2777777777778, 21.0526315789474, 27.3239436619718, 22.1105527638191, 22.9268292682927, 24.1071428571429, 29.5522388059701, 25.6267409470752, 31.8611987381703, 30.1994301994302, 26.2337662337662, 23.9608801955990, 31.6614420062696, 24.8730964467005, 17.8010471204188, 24.2506811989101, 36.1445783132530, 18.8888888888889, 25.4847645429363, 24.9343832020997, 27.1767810026385, 24.3176178660050, 22.7528089887640, 18.4696569920844, 23.7430167597765, 25.2873563218391, 22.6256983240223, 19.1601049868766, 21.4470284237726, 28.2321899736148, 23.0769230769231], [13.4730538922156, 11.8733509234828, 9.86301369863014, 9.62566844919786, 10.8433734939759, 13.9896373056995, 15.5172413793103, 12.0967741935484, 12.5000000000000, 11.2781954887218, 10.1408450704225, 9.04522613065327, 10.9756097560976, 13.3928571428571, 10.7462686567164, 12.5348189415042, 11.3564668769716, 10.2564102564103, 9.35064935064935, 11.0024449877751, 11.2852664576803, 11.4213197969543, 14.1361256544503, 12.2615803814714, 14.4578313253012, 12.5000000000000, 12.4653739612188, 9.44881889763779, 9.49868073878628, 8.93300248138958, 12.6404494382022, 11.8733509234828, 12.5698324022346, 12.9310344827586, 12.5698324022346, 14.1732283464567, 11.6279069767442, 9.49868073878628, 12.3626373626374]]

def seq_length(P):
    for k in P._structures:
        seq_length = len(P._structures[k])
        break
    return seq_length
def get_all_polys(sobj_names):
    Polytopes = {}
    for sname in sobj_names:
        Polytopes[sname] = RNAPolytope.construct_from_file(sname)
    return Polytopes

print timestamp()
n_facets = defaultdict(list)
for sname in Polytopes:
    length = seq_length(Polytopes[sname])
    n_facets[length].append(len(Polytopes[sname].faces(3)))
print timestamp()
# 17 minutes

[2017-08-10 19:21:25.056500] [2017-08-10 19:38:38.601972]

n_facets

defaultdict(<type 'list'>, {72: [2106, 2177, 1696, 2036], 73: [2247, 1826, 1848, 1917, 2026, 1752, 1644, 1890, 2148, 1806, 2212, 1930], 74: [1780, 2139, 2245, 2161, 1957, 2019, 2100, 2130], 75: [2201, 2281], 76: [2121, 2068, 2084, 2176, 2347, 1985], 77: [2529, 2145, 2123, 2506, 2343, 2382], 79: [2058]})

n_facets

defaultdict(<type 'list'>, {76: [2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176, 2176]})

for k in P._structures:
    seq_length = len(P._structures[k])
print seq_length

76

vx_totals={}
print timestamp()
for sname in sobj_names:
    P=RNAPolytope.construct_from_file(sname)
    vx_totals[sname] = len(P.vertices_list())
print timestamp()

[2017-08-10 18:41:52.651620] [2017-08-10 18:44:50.942809]

pct_d1 = []
lost = []
for sname in sobj_names:
    pct_d1.append(100.0*float(totals[sname])/vx_totals[sname])
    lost.append(vx_totals[sname]-totals[sname])

print pct_d1
print avg(pct_d1)
print std(pct_d1)
print med(pct_d1)
print lost
print avg(lost)
print std(lost)
print med(lost)

[70.0819672131148, 72.2222222222222, 71.1051930758988, 74.0786240786241, 75.0000000000000, 75.9112519809826, 78.0000000000000, 75.8379888268156, 69.6132596685083, 74.6206896551724, 66.3143989431968, 70.1680672268908, 73.3236151603499, 77.5474956822107, 71.4677640603567, 79.5562599049128, 74.8129675810474, 71.0937500000000, 70.2341137123746, 76.2032085561497, 73.7878787878788, 73.3516483516484, 72.5903614457831, 76.5536723163842, 70.2325581395349, 72.0990873533246, 77.7927321668910, 74.6121297602257, 73.1323722149410, 67.9566563467492, 78.9687924016282, 73.9946380697051, 75.2659574468085, 67.3352435530086, 77.3121387283237, 71.6772151898734, 75.1085383502171, 68.4082624544350, 74.2088607594937] 73.3738867022 3.22618211829 73.7878787879 [219, 190, 217, 211, 173, 152, 143, 173, 220, 184, 255, 213, 183, 130, 208, 129, 202, 185, 178, 178, 173, 194, 182, 166, 192, 214, 165, 180, 205, 207, 155, 194, 186, 228, 157, 179, 172, 260, 163] 187.564102564 28.8584475228 184.0

# Average percent of bounded regions
pct_d1 = [70.0819672131148, 72.2222222222222, 71.1051930758988, 74.0786240786241, 75.0000000000000, 75.9112519809826, 78.0000000000000, 75.8379888268156, 69.6132596685083, 74.6206896551724, 66.3143989431968, 70.1680672268908, 73.3236151603499, 77.5474956822107, 71.4677640603567, 79.5562599049128, 74.8129675810474, 71.0937500000000, 70.2341137123746, 76.2032085561497, 73.7878787878788, 73.3516483516484, 72.5903614457831, 76.5536723163842, 70.2325581395349, 72.0990873533246, 77.7927321668910, 74.6121297602257, 73.1323722149410, 67.9566563467492, 78.9687924016282, 73.9946380697051, 75.2659574468085, 67.3352435530086, 77.3121387283237, 71.6772151898734, 75.1085383502171, 68.4082624544350, 74.2088607594937]

lost

[219, 190, 217, 211, 173, 152, 143, 173, 220, 184, 255, 213, 183, 130, 208, 129, 202, 185, 178, 178, 173, 194, 182, 166, 192, 214, 165, 180, 205, 207, 155, 194, 186, 228, 157, 179, 172, 260, 163]

print timestamp()
boundeds = {}
for sname in sobj_names:
    boundeds[sname] = get_bounded_partitions(sname)
print timestamp()

[2017-08-10 18:36:04.920589] [2017-08-10 18:36:34.805950] Unbounded only in ab: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 0, 22, 177/5); (351/20, 29/15, -509/40) [2017-08-10 18:36:34.806913] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 12, 22, 189/5); (67138/1875, -8047/3750, -41126/1875) [2017-08-10 18:36:34.808069] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 8, 22, 169/5); (1637/72, -947/1800, -26791/1800) [2017-08-10 18:36:34.809573] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 10, 22, 71/2); (3837/140, -129/140, -9369/560) [2017-08-10 18:36:34.810384] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 3, 22, 167/5); (1801/80, -1/40, -1197/80) [2017-08-10 18:36:34.811937] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 2, 22, 167/5); (4381/220, 113/220, -12213/880) [2017-08-10 18:36:34.813137] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 1, 22, 171/5); (2491/135, 269/270, -3623/270) [2017-08-10 18:36:34.815143] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 11, 22, 73/2); (42277/1560, -67/60, -51457/3120) [2017-08-10 18:36:34.816170] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 7, 22, 168/5); (7957/360, -5/36, -5189/360) [2017-08-10 18:37:17.349598] Unbounded only in ab: sobjs/nocardia.farcinica.ifm.10152_AP006618_az.sobj; A vertex at (8, 0, 25, 264/5); (11423/800, 37/8, -1169/100) [2017-08-10 18:38:07.936114] Unbounded only in ab: sobjs/hyperthermus.butylicus.dsm.5456_CP000493_ad.sobj; A vertex at (8, 0, 25, 4691/100); (1203/100, 1433/140, -67/5) [2017-08-10 18:38:15.951254] Unbounded only in ab: sobjs/uncultured.gamma.proteobacterium_AY934387_a.sobj; A vertex at (8, 0, 25, 549/10); (1009/80, 511/80, -175/16) [2017-08-10 18:38:31.275948] Unbounded only in ab: sobjs/sulfolobus.tokodaii.str..7_BA000023_ak.sobj; A vertex at (8, 0, 25, 49); (13421/600, 1817/600, -2871/200) [2017-08-10 18:38:38.864601]

for n in boundeds:
    tot = totals

%md

There are only **8 regions that are strips** on a and all of them occur at **synechococcus**.

```
[2017-08-10 01:34:41.059542] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 12, 22, 189/5); (67138/1875, -8047/3750, -41126/1875)
[2017-08-10 01:34:41.060729] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 8, 22, 169/5); (1637/72, -947/1800, -26791/1800)
[2017-08-10 01:34:41.062320] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 10, 22, 71/2); (3837/140, -129/140, -9369/560)
[2017-08-10 01:34:41.063201] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 3, 22, 167/5); (1801/80, -1/40, -1197/80)
[2017-08-10 01:34:41.064876] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 2, 22, 167/5); (4381/220, 113/220, -12213/880)
[2017-08-10 01:34:41.065926] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 1, 22, 171/5); (2491/135, 269/270, -3623/270)
[2017-08-10 01:34:41.067180] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 11, 22, 73/2); (42277/1560, -67/60, -51457/3120)
[2017-08-10 01:34:41.068276] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 7, 22, 168/5); (7957/360, -5/36, -5189/360)
```

There are only **5 regions that are unbounded on ab** and there is one on each of the following sequences.

```
[2017-08-10 01:34:41.058388] Unbounded only in ab: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 0, 22, 177/5); (351/20, 29/15, -509/40)
[2017-08-10 01:35:21.704100] Unbounded only in ab: sobjs/nocardia.farcinica.ifm.10152_AP006618_az.sobj; A vertex at (8, 0, 25, 264/5); (11423/800, 37/8, -1169/100)
[2017-08-10 01:36:18.094037] Unbounded only in ab: sobjs/hyperthermus.butylicus.dsm.5456_CP000493_ad.sobj; A vertex at (8, 0, 25, 4691/100); (1203/100, 1433/140, -67/5)
[2017-08-10 01:36:26.910393] Unbounded only in ab: sobjs/uncultured.gamma.proteobacterium_AY934387_a.sobj; A vertex at (8, 0, 25, 549/10); (1009/80, 511/80, -175/16)
[2017-08-10 01:36:42.061085] Unbounded only in ab: sobjs/sulfolobus.tokodaii.str..7_BA000023_ak.sobj; A vertex at (8, 0, 25, 49); (13421/600, 1817/600, -2871/200)
```

There are only 8 regions that are strips on a and all of them occur at synechococcus.

[2017-08-10 01:34:41.059542] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 12, 22, 189/5); (67138/1875, -8047/3750, -41126/1875)
[2017-08-10 01:34:41.060729] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 8, 22, 169/5); (1637/72, -947/1800, -26791/1800)
[2017-08-10 01:34:41.062320] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 10, 22, 71/2); (3837/140, -129/140, -9369/560)
[2017-08-10 01:34:41.063201] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 3, 22, 167/5); (1801/80, -1/40, -1197/80)
[2017-08-10 01:34:41.064876] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 2, 22, 167/5); (4381/220, 113/220, -12213/880)
[2017-08-10 01:34:41.065926] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 1, 22, 171/5); (2491/135, 269/270, -3623/270)
[2017-08-10 01:34:41.067180] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 11, 22, 73/2); (42277/1560, -67/60, -51457/3120)
[2017-08-10 01:34:41.068276] Unbounded only in a: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 7, 22, 168/5); (7957/360, -5/36, -5189/360)

There are only 5 regions that are unbounded on ab and there is one on each of the following sequences.

[2017-08-10 01:34:41.058388] Unbounded only in ab: sobjs/synechococcus.sp..wh.8102_BX569691.sobj; A vertex at (7, 0, 22, 177/5); (351/20, 29/15, -509/40)
[2017-08-10 01:35:21.704100] Unbounded only in ab: sobjs/nocardia.farcinica.ifm.10152_AP006618_az.sobj; A vertex at (8, 0, 25, 264/5); (11423/800, 37/8, -1169/100)
[2017-08-10 01:36:18.094037] Unbounded only in ab: sobjs/hyperthermus.butylicus.dsm.5456_CP000493_ad.sobj; A vertex at (8, 0, 25, 4691/100); (1203/100, 1433/140, -67/5)
[2017-08-10 01:36:26.910393] Unbounded only in ab: sobjs/uncultured.gamma.proteobacterium_AY934387_a.sobj; A vertex at (8, 0, 25, 549/10); (1009/80, 511/80, -175/16)
[2017-08-10 01:36:42.061085] Unbounded only in ab: sobjs/sulfolobus.tokodaii.str..7_BA000023_ak.sobj; A vertex at (8, 0, 25, 49); (13421/600, 1817/600, -2871/200)

totals=defaultdict(lambda:0)
for n in boundeds:
    totals[n] = sum([boundeds[n][k] for k in boundeds[n] if len(k)<4])

boundeds_pcts = {}
for n in boundeds:
    boundeds_pcts[n]={}
    for k in boundeds[n]:
        if len(k)>4:
            boundeds_pcts[n][k] = boundeds[n][k]
            continue
        boundeds_pcts[n][k] = (100.0*boundeds[n][k])/totals[n]

for i,n in enumerate(boundeds_pcts):
    print n[:20],
    for k in ['bounded','unbounded','ac','abc','b','ab','a']:
        if k not in boundeds[n]:
            continue
        print "{0}: {1: 6.3f};".format(k,float(boundeds_pcts[n][k])),
    print
    #boundeds_pcts[n]

sobjs/chlorobium.pha bounded: 60.220; unbounded: 39.780; ac: 23.297; abc: 12.747; b: 3.736; sobjs/homo.sapiens_X bounded: 64.870; unbounded: 35.130; ac: 17.964; abc: 12.974; b: 4.192; sobjs/solibacter.usi bounded: 65.047; unbounded: 34.953; ac: 20.561; abc: 11.776; b: 2.617; sobjs/uncultured.gam bounded: 60.777; unbounded: 39.223; ac: 24.558; abc: 10.954; b: 3.534; ab: 0.177; sobjs/pyrococcus.aby bounded: 66.842; unbounded: 33.158; ac: 17.719; abc: 11.404; b: 4.035; sobjs/yersinia.pesti bounded: 63.352; unbounded: 36.648; ac: 21.179; abc: 11.602; b: 3.867; sobjs/yarrowia.lipol bounded: 61.276; unbounded: 38.724; ac: 20.729; abc: 13.212; b: 4.784; sobjs/sulfolobus.tok bounded: 59.345; unbounded: 40.655; ac: 23.699; abc: 13.487; b: 3.276; ab: 0.193; sobjs/euglena.longa_ bounded: 59.138; unbounded: 40.862; ac: 24.846; abc: 11.499; b: 4.517; sobjs/bordetella.par bounded: 65.217; unbounded: 34.783; ac: 18.715; abc: 11.909; b: 4.159; sobjs/nitrosomonas.e bounded: 60.234; unbounded: 39.766; ac: 24.172; abc: 11.891; b: 3.704; sobjs/rickettsia.con bounded: 62.979; unbounded: 37.021; ac: 17.872; abc: 13.404; b: 5.745; sobjs/uncultured.alp bounded: 65.174; unbounded: 34.826; ac: 20.398; abc: 11.277; b: 3.151; sobjs/methanospirill bounded: 60.044; unbounded: 39.956; ac: 23.841; abc: 12.141; b: 3.974; sobjs/wigglesworthia bounded: 59.382; unbounded: 40.618; ac: 24.945; abc: 13.466; b: 2.208; sobjs/arabidopsis.th bounded: 61.233; unbounded: 38.767; ac: 23.260; abc: 12.127; b: 3.380; sobjs/synechococcus. bounded: 55.621; unbounded: 44.379; ac: 28.994; abc: 10.256; b: 3.353; ab: 0.197; a: 1.578; sobjs/uncultured.alt bounded: 60.701; unbounded: 39.299; ac: 24.908; abc: 11.624; b: 2.768; sobjs/magnetococcus. bounded: 62.903; unbounded: 37.097; ac: 22.043; abc: 11.649; b: 3.405; sobjs/pelobacter.car bounded: 65.052; unbounded: 34.948; ac: 19.896; abc: 11.938; b: 3.114; sobjs/nephroselmis.o bounded: 56.290; unbounded: 43.710; ac: 28.785; abc: 11.727; b: 3.198; sobjs/tropheryma.whi bounded: 60.501; unbounded: 39.499; ac: 23.121; abc: 11.946; b: 4.432; sobjs/marchantia.pol bounded: 64.981; unbounded: 35.019; ac: 19.476; abc: 11.236; b: 4.307; sobjs/homo.sapiens_A bounded: 60.145; unbounded: 39.855; ac: 25.181; abc: 10.688; b: 3.986; sobjs/scenedesmus.ob bounded: 51.429; unbounded: 48.571; ac: 35.238; abc: 11.190; b: 2.143; sobjs/uncultured.gam bounded: 64.542; unbounded: 35.458; ac: 20.319; abc: 11.952; b: 3.187; sobjs/lactobacillus. bounded: 59.109; unbounded: 40.891; ac: 24.899; abc: 11.741; b: 4.251; sobjs/drosophila.mel bounded: 65.099; unbounded: 34.901; ac: 20.615; abc: 11.392; b: 2.893; sobjs/carboxydotherm bounded: 60.259; unbounded: 39.741; ac: 23.475; abc: 12.569; b: 3.697; sobjs/uncultured.mar bounded: 64.831; unbounded: 35.169; ac: 19.716; abc: 12.078; b: 3.375; sobjs/uncultured.alp bounded: 61.310; unbounded: 38.690; ac: 23.214; abc: 11.508; b: 3.968; sobjs/chlamydophila. bounded: 63.920; unbounded: 36.080; ac: 17.595; abc: 15.367; b: 3.118; sobjs/haloferax.volc bounded: 60.558; unbounded: 39.442; ac: 23.506; abc: 12.550; b: 3.386; sobjs/thermoplasma.v bounded: 62.188; unbounded: 37.812; ac: 23.800; abc: 10.749; b: 3.263; sobjs/uncultured.thi bounded: 61.795; unbounded: 38.205; ac: 22.338; abc: 13.152; b: 2.714; sobjs/methanosarcina bounded: 65.543; unbounded: 34.457; ac: 18.539; abc: 11.423; b: 4.494; sobjs/nocardia.farci bounded: 66.167; unbounded: 33.833; ac: 20.167; abc: 10.833; b: 2.667; ab: 0.167; sobjs/hyperthermus.b bounded: 60.997; unbounded: 39.003; ac: 24.914; abc: 11.340; b: 2.577; ab: 0.172; sobjs/lactobacillus. bounded: 59.336; unbounded: 40.664; ac: 23.237; abc: 13.071; b: 4.357;

out_bound = defaultdict(list)
for n in boundeds_pcts:
    for k in boundeds_pcts[n]:
        out_bound[k].append(boundeds_pcts[n][k])

values = []
labels = ['bounded','unbounded','ac','abc','b']
for k in labels:
    values.append([float(x) for x in out_bound[k]])
generic_boxplot(values,labels)

:)

Error in lines 5-5 Traceback (most recent call last): File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/smc_sagews/sage_server.py", line 995, in execute exec compile(block+'\n', '', 'single') in namespace, locals File "", line 1, in <module> File "", line 13, in generic_boxplot File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/matplotlib/pyplot.py", line 2710, in boxplot manage_xticks=manage_xticks, data=data) File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/matplotlib/__init__.py", line 1812, in inner return func(ax, *args, **kwargs) File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/matplotlib/axes/_axes.py", line 3323, in boxplot manage_xticks=manage_xticks) File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/matplotlib/axes/_axes.py", line 3598, in bxp elif len(widths) != N: TypeError: object of type 'sage.rings.real_mpfr.RealLiteral' has no len()

values

[[60.21978021978022, 64.87025948103792, 65.04672897196262, 60.7773851590106, 66.84210526315789, 63.35174953959484, 61.275626423690206, 59.34489402697495, 59.13757700205339, 65.21739130434783, 60.23391812865497, 62.97872340425532, 65.17412935323384, 60.04415011037528, 59.381898454746135, 61.232604373757454, 55.62130177514793, 60.701107011070114, 62.903225806451616, 65.05190311418686, 56.28997867803838, 60.5009633911368, 64.9812734082397, 60.14492753623188, 51.42857142857143, 64.54183266932272, 59.10931174089069, 65.0994575045208, 60.258780036968574, 64.83126110124334, 61.30952380952381, 63.919821826280625, 60.55776892430279, 62.18809980806142, 61.79540709812109, 65.54307116104869, 66.16666666666667, 60.99656357388316, 59.33609958506224], [39.78021978021978, 35.12974051896207, 34.953271028037385, 39.2226148409894, 33.1578947368421, 36.64825046040516, 38.724373576309794, 40.65510597302505, 40.86242299794661, 34.78260869565217, 39.76608187134503, 37.02127659574468, 34.82587064676617, 39.95584988962472, 40.618101545253865, 38.767395626242546, 44.37869822485207, 39.298892988929886, 37.096774193548384, 34.94809688581315, 43.71002132196162, 39.4990366088632, 35.0187265917603, 39.85507246376812, 48.57142857142857, 35.45816733067729, 40.89068825910931, 34.90054249547921, 39.741219963031426, 35.16873889875666, 38.69047619047619, 36.080178173719375, 39.44223107569721, 37.81190019193858, 38.20459290187891, 34.45692883895131, 33.833333333333336, 39.00343642611684, 40.66390041493776], [23.296703296703296, 17.964071856287426, 20.560747663551403, 24.558303886925795, 17.719298245614034, 21.17863720073665, 20.72892938496583, 23.69942196531792, 24.84599589322382, 18.714555765595463, 24.171539961013647, 17.872340425531913, 20.398009950248756, 23.841059602649008, 24.944812362030905, 23.26043737574553, 28.994082840236686, 24.907749077490774, 22.043010752688172, 19.896193771626297, 28.784648187633262, 23.121387283236995, 19.475655430711612, 25.181159420289855, 35.23809523809524, 20.318725099601593, 24.898785425101213, 20.61482820976492, 23.475046210720887, 19.715808170515096, 23.214285714285715, 17.59465478841871, 23.50597609561753, 23.80038387715931, 22.338204592901878, 18.53932584269663, 20.166666666666668, 24.914089347079038, 23.236514522821576], [12.747252747252746, 12.974051896207586, 11.77570093457944, 10.954063604240282, 11.403508771929825, 11.602209944751381, 13.211845102505695, 13.48747591522158, 11.498973305954825, 11.909262759924385, 11.89083820662768, 13.404255319148936, 11.276948590381426, 12.141280353200884, 13.465783664459162, 12.127236580516898, 10.256410256410257, 11.623616236162361, 11.64874551971326, 11.937716262975778, 11.727078891257996, 11.946050096339114, 11.235955056179776, 10.68840579710145, 11.19047619047619, 11.952191235059761, 11.740890688259109, 11.39240506329114, 12.569316081330868, 12.078152753108348, 11.507936507936508, 15.367483296213809, 12.549800796812749, 10.748560460652591, 13.152400835073069, 11.423220973782772, 10.833333333333334, 11.34020618556701, 13.070539419087137], [3.7362637362637363, 4.191616766467066, 2.616822429906542, 3.5335689045936394, 4.035087719298246, 3.867403314917127, 4.783599088838269, 3.275529865125241, 4.517453798767967, 4.158790170132325, 3.7037037037037037, 5.74468085106383, 3.1509121061359866, 3.9735099337748343, 2.207505518763797, 3.3797216699801194, 3.353057199211045, 2.7675276752767526, 3.4050179211469533, 3.114186851211073, 3.1982942430703623, 4.431599229287091, 4.307116104868914, 3.9855072463768115, 2.142857142857143, 3.187250996015936, 4.251012145748988, 2.8933092224231465, 3.6968576709796674, 3.374777975133215, 3.9682539682539684, 3.11804008908686, 3.3864541832669324, 3.2629558541266794, 2.7139874739039667, 4.49438202247191, 2.6666666666666665, 2.577319587628866, 4.356846473029045]]

set([tuple(sorted(boundeds[n].keys())) for n in boundeds])

set([('', 'ab', 'abc', 'ac', 'b', 'bounded', 'unbounded'), ('', 'abc', 'ac', 'b', 'bounded', 'unbounded'), ('', 'a', 'ab', 'abc', 'ac', 'b', 'bounded', 'unbounded')])

merge_bound = defaultdict(lambda:0)
for n in boundeds:
    for k in boundeds[n]:
        merge_bound[k]+= boundeds[n][k]

print merge_bound

defaultdict(<function <lambda> at 0x7f1de7a427d0>, {'': 12459, 'a': 8, 'ac': 4529, 'b': 718, 'unbounded': 1491.59416112840, 'abc': 2405, 'bounded': 2408.40583887160, 'ab': 5})

sl.dimension()

3

len(P.vertices_list())
len(P.d1_slices())

732 513

get_shade = lambda x,M,m: (x-m)/(M-m)

def get_all_hwd_tuples(sobj_names,by_val='x'):
    hwd_tuples = []

    if by_val=='x':
        ind = 0
    elif by_val=='y':
        ind = 1
    elif by_val=='z':
        ind = 2
    else:
        raise ValueError("Unknown by_val, must be one of 'x', 'y', or 'z'.")

    polytopes = [RNAPolytope.construct_from_file(sname) for sname in sobj_names]

    max_val = max([sl.original_vertex[ind] for P in polytopes for sl in P.d1_slices() if is_bounded(sl)])
    min_val = min([sl.original_vertex[ind] for P in polytopes for sl in P.d1_slices() if is_bounded(sl)])

    hwd_tuples = {}
    for P in polytopes:
        for sl in P.d1_slices():
            if is_bounded(sl):
                val = sl.original_vertex[ind]
                hwd_tuple = get_slice_width_height_depth(sl)
                if val in hwd_tuples:
                    hwd_tuples[val].append(hwd_tuple)
                else:
                    hwd_tuples[val]=[hwd_tuple]
    new_hwd_tuples = {}
    for val in hwd_tuples:
        try:
            shade = get_shade(val,max_val,min_val)
        except ZeroDivisionError as e:
            print max_val,min_val
            print sorted(hwd_tuples.keys())
            raise e
        new_hwd_tuples[shade] = hwd_tuples[val]
    return new_hwd_tuples


def get_hwd_tuples(P,by_val='x'):
    hwd_tuples = []

    if by_val=='x':
        ind = 0
    elif by_val=='y':
        ind = 1
    elif by_val=='z':
        ind = 2
    else:
        raise ValueError("Unknown by_val, must be one of 'x', 'y', or 'z'.")
    max_val = max([sl.original_vertex[ind] for sl in P.d1_slices() if is_bounded(sl)])
    min_val = min([sl.original_vertex[ind] for sl in P.d1_slices() if is_bounded(sl)])

    hwd_tuples = {}

    for sl in P.d1_slices():
        if is_bounded(sl):
            val = sl.original_vertex[ind]
            hwd_tuple = get_slice_width_height_depth(sl)
            if val in hwd_tuples:
                hwd_tuples[val].append(hwd_tuple)
            else:
                hwd_tuples[val]=[hwd_tuple]
    new_hwd_tuples = {}
    for val in hwd_tuples:
        try:
            shade = get_shade(val,max_val,min_val)
        except ZeroDivisionError as e:
            print max_val,min_val
            print sorted(hwd_tuples.keys())
            raise e
        new_hwd_tuples[shade] = hwd_tuples[val]
    return new_hwd_tuples


def plot_hwds(hwds,start=0):
    pl = Graphics()
    keys = sorted(hwds.keys())
    for shade in keys[start:]:
        pl += point3d(hwds[shade], color=my_hot(shade))
    return pl

# 5 minutes
# def get_all_hwd_tuples():
#     hwd_tuples = []
#     n=len(sobj_names)
#
#     for i,sname in enumerate(sobj_names):
#         print timestamp('[{0: 3d}/{1: 3d}] Started {2}'.format(i+1,n,sname))
#         P = RNAPolytope.construct_from_file(sname)
#         hwd_tuples.extend(get_hwd_tuples(P))
#         print timestamp('[{0: 3d}/{1: 3d}] Finishd {2}'.format(i+1,n,sname))
#     return hwd_tuples

hwds = get_all_hwd_tuples(sobj_names)

flat_hwd = []
for k in hwds:
    flat_hwd.extend(hwds[k])

flat1_hwd = []
for k in sorted(hwds.keys()):
    if k>0:
        flat1_hwd.extend([(x,z) for x,y,z in hwds[k]])
    print "{0}: {1}".format(k,len(hwds[k]))
flat1_float = [(float(x),float(y)) for x,y in flat1_hwd]

0: 3608 1/7: 2783 2/7: 2187 3/7: 1570 4/7: 1154 5/7: 733 6/7: 347 1: 77

import pickle
def save_obj(obj, name):
    with open(name, 'wb') as f:
        pickle.dump(obj, f, pickle.HIGHEST_PROTOCOL)

def load_obj(name):
    with open(name, 'rb') as f:
        return pickle.load(f)
save_obj(flat1_float,'ac_ratios.pkl')

cvar = [c for a,b,c in flat_hwd]
float(min(cvar))
float(max(cvar))
avg(cvar)
std(cvar)
float(median(cvar))

0.01818181818181818 30.9 2.6067763501299224 2.6729232154816196 1.8333333333333333

avar = [a for a,b,c in flat_hwd]
float(min(avar))
float(max(avar))
avg(avar)
std(avar)
float(median(avar))

0.03 498.25 28.259385014501277 50.445283479276704 5.8

bvar = [b for a,b,c in flat_hwd]
float(min(bvar))
float(max(bvar))
avg(bvar)
std(bvar)
float(median(bvar))

0.0011111111111111111 6.1 0.44068681304337654 0.51112271239128726 0.3

min([y for x,y,z in hwds])
max([y for x,y,z in hwds])
float(min([y for x,y,z in hwds]))
float(max([y for x,y,z in hwds]))

1/900 61/10 0.0011111111111111111 6.1

#point3d(hwds)

#Q=RNAPolytope.construct_from_file(sobj_names[0])
Qhwds= get_hwd_tuples(Q)
#point3d(hwds)

Xhwds=get_all_hwd_tuples(sobj_names[:3],'x')

#plot_hwds(Rhwds,start=1)

# ~5 min (uncommented version)
t = timestamp()
print t
# Xhwds=get_all_hwd_tuples(sobj_names[:10],'x')
# k =  len(Xhwds.keys())
# print "{0} different values: {1}".format(k,sorted(Xhwds.keys()))
# plot_hwds(Xhwds,start=0)
# plot_hwds(Xhwds,start=1)
# plot_hwds(Xhwds,start=int(k/5))
# Yhwds=get_all_hwd_tuples(sobj_names[:10],'y')
# k =  len(Yhwds.keys())
# print "{0} different values: {1}".format(k,sorted(Yhwds.keys()))
# plot_hwds(Yhwds,start=0)
# plot_hwds(Yhwds,start=1)
# plot_hwds(Yhwds,start=int(k/4))
# plot_hwds(Yhwds,start=int(k/3))
# plot_hwds(Yhwds,start=int(k/2))
# Zhwds=get_all_hwd_tuples(sobj_names[:10],'z')
# k =  len(Zhwds.keys())
# print "{0} different values: {1}".format(k,sorted(Zhwds.keys()))
# plot_hwds(Zhwds,start=0)
# plot_hwds(Zhwds,start=1)
# plot_hwds(Zhwds,start=int(k/4))
# plot_hwds(Zhwds,start=int(k/3))
# plot_hwds(Zhwds,start=int(k/2))
print timestamp()

[2017-08-20 23:03:12.990176] [2017-08-20 23:03:12.991073]

plot_hwds(Rhwds,start=1)

Qhwds.keys()

[0, 1, 1/2, 1/3, 2/3, 1/6, 5/6]

for sname in sobj_names[:3]:
    R = RNAPolytope.construct_from_file(sname)
    Rhwds = get_hwd_tuples(R)
    print sname
    plot_hwds(Rhwds,start=1)

get_slice_width_height_depth(sl)

[4283/30, 19/30, 11/3]

get_slice_width_height_depth_intervals(sl)

[(787/30, 169), (-12/5, -53/30), (-169/10, -397/30)]

min([x for x,y,z in sl.vertices()])
max([x for x,y,z in sl.vertices()])
min([y for x,y,z in sl.vertices()])
max([y for x,y,z in sl.vertices()])
min([z for x,y,z in sl.vertices()])
max([z for x,y,z in sl.vertices()])

787/30 169 -12/5 -53/30 -169/10 -397/30

plot_slice(sl)

accs,structs,signs=compute_subopt_similarities(names[0])

[2017-08-02 05:54:43.469318] nitrosomonas.europaea.atcc.19718_BX321863_b: max_subopts=5, avg_subopts=1.23781676413; std_subopts=0.6228062549

stat_accs = []
minsigs = []
minstructs = []
for name in names:
    #name = 'wigglesworthia.glossinidia.endosymbiont.of.glossina.brevipalpis_BA000021_ad'
    accs,structs,signs=compute_subopt_similarities(name)
    minacc = min(accs.values())
    avgacc = avg(accs.values())
    stdacc = std(accs.values())
    medacc = std(accs.values())
    stat_accs.append((minacc,avgacc,stdacc,medacc))

    minsig = [s for s in accs if accs[s]==minacc]
    minsigs.append(minsig)
    min_seq_struct = []
    for sig in minsig:
        min_sig_struct=[]
        for s,t in combinations(structs[sig],2):
            if f_measure(s,t)==minacc:
                min_sig_struct.append((s,t))
            elif f_measure(t,s)==minacc:
                min_sig_struct.append((t,s))
        min_seq_struct.append(min_sig_struct)
    minstructs.append(min_seq_struct)
    ss=structs[minsig[0]]
    print name
    print minacc
    print minsig,len(ss)
    print avg(accs.values())
    print std(accs.values())
    print "*"*50
    #break

[2017-08-02 15:23:35.256941] nitrosomonas.europaea.atcc.19718_BX321863_b: max_subopts=5, avg_subopts=1.23781676413; std_subopts=0.6228062549 nitrosomonas.europaea.atcc.19718_BX321863_b 0.538461538462 [(6, 25, 18, (133, 10))] 2 0.834658632373 0.107724018943 ************************************************** [2017-08-02 15:23:35.285255] lactobacillus.casei.atcc.334_CP000423_j: max_subopts=4, avg_subopts=1.10728744939; std_subopts=0.358002118099 lactobacillus.casei.atcc.334_CP000423_j 0.228571428571 [(7, 8, 21, (47, 2))] 3 0.818226996303 0.155252521689 ************************************************** [2017-08-02 15:23:35.355989] methanosarcina.barkeri.str..fusaro_CP000099_bi: max_subopts=15, avg_subopts=1.37827715356; std_subopts=1.13967274278 methanosarcina.barkeri.str..fusaro_CP000099_bi 0.235294117647 [(6, 6, 19, (35, 2))] 2 0.810933070304 0.149287342991 ************************************************** [2017-08-02 15:23:35.437880] uncultured.alpha.proteobacterium_AY933788: max_subopts=12, avg_subopts=1.42122719735; std_subopts=1.00166386608 uncultured.alpha.proteobacterium_AY933788 0.272727272727 [(5, 37, 15, (187, 10))] 2 0.804289185104 0.142093535926 ************************************************** [2017-08-02 15:23:35.522373] tropheryma.whipplei.str..twist_AE014184_h: max_subopts=8, avg_subopts=1.49132947977; std_subopts=0.955373185828 tropheryma.whipplei.str..twist_AE014184_h 0.291666666667 [(4, 3, 12, (-111, 5))] 2 0.873287591355 0.111832176641 ************************************************** [2017-08-02 15:23:35.584935] uncultured.thiotrichales.bacterium_AY934254_a: max_subopts=7, avg_subopts=1.27557411273; std_subopts=0.64232130098 uncultured.thiotrichales.bacterium_AY934254_a 0.461538461538 [(5, 22, 16, (77, 10))] 2 0.863202920273 0.122686172757 ************************************************** [2017-08-02 15:23:35.673408] synechococcus.sp..wh.8102_BX569691: max_subopts=12, avg_subopts=1.57790927022; std_subopts=1.25274655561 synechococcus.sp..wh.8102_BX569691 0.181818181818 [(7, 12, 21, (243, 10))] 4 0.843787946158 0.0999454149297 ************************************************** [2017-08-02 15:23:35.740494] yersinia.pestis.kim_AE009952_ae: max_subopts=5, avg_subopts=1.12154696133; std_subopts=0.407068983779 yersinia.pestis.kim_AE009952_ae 0.1875 [(4, 16, 12, (-129, 10))] 3 0.768702893294 0.192479942686 ************************************************** [2017-08-02 15:23:35.793947] uncultured.alpha.proteobacterium_DQ093144_a: max_subopts=6, avg_subopts=1.27380952381; std_subopts=0.729500652374 uncultured.alpha.proteobacterium_DQ093144_a 0.325581395349 [(1, 0, 8, (5, 2))] 2 0.789794159242 0.165393594841 ************************************************** [2017-08-02 15:23:35.836995] carboxydothermus.hydrogenoformans.z-2901_CP000141_af: max_subopts=5, avg_subopts=1.11829944547; std_subopts=0.413337785292 carboxydothermus.hydrogenoformans.z-2901_CP000141_af 0.176470588235 [(8, 13, 24, (167, 5))] 5 0.805549974678 0.169076757089 ************************************************** [2017-08-02 15:23:35.872978] haloferax.volcanii_HALTRK1: max_subopts=4, avg_subopts=1.14342629482; std_subopts=0.408264490214 haloferax.volcanii_HALTRK1 0.162162162162 [(8, 3, 25, (493, 10))] 2 0.763001477859 0.201427585696 ************************************************** [2017-08-02 15:23:35.917407] homo.sapiens_X04779: max_subopts=5, avg_subopts=1.22155688623; std_subopts=0.506262548212 homo.sapiens_X04779 0.242424242424 [(3, 0, 15, (39, 1))] 3 0.840567749448 0.151216389289 ************************************************** [2017-08-02 15:23:35.966903] arabidopsis.thaliana_AC006340: max_subopts=6, avg_subopts=1.29423459245; std_subopts=0.688521740769 arabidopsis.thaliana_AC006340 0.358974358974 [(6, 2, 18, (13, 2)), (6, 3, 18, (51, 10))] 2 0.837990094411 0.138631934747 ************************************************** [2017-08-02 15:23:36.004735] chlamydophila.pneumoniae.tw-183_AE017159_j: max_subopts=9, avg_subopts=1.25389755011; std_subopts=0.665831160483 chlamydophila.pneumoniae.tw-183_AE017159_j 0.205128205128 [(7, 0, 21, (63, 2))] 2 0.845632728425 0.13522223029 ************************************************** [2017-08-02 15:23:36.068984] thermoplasma.volcanium.gss1_BA000011_b: max_subopts=11, avg_subopts=1.25719769674; std_subopts=0.742894051181 thermoplasma.volcanium.gss1_BA000011_b 0.28 [(2, 4, 6, (-351, 10))] 2 0.834635603067 0.142729785488 ************************************************** [2017-08-02 15:23:36.113078] uncultured.gamma.proteobacterium_AY934184_a: max_subopts=3, avg_subopts=1.12151394422; std_subopts=0.39817224196 uncultured.gamma.proteobacterium_AY934184_a 0.444444444444 [(4, 15, 16, (331, 10))] 2 0.844670013225 0.128834330911 ************************************************** [2017-08-02 15:23:36.160945] nocardia.farcinica.ifm.10152_AP006618_az: max_subopts=6, avg_subopts=1.165; std_subopts=0.494410423299 nocardia.farcinica.ifm.10152_AP006618_az 0.216216216216 [(6, 9, 18, (32, 5))] 2 0.785856379183 0.177105262683 ************************************************** [2017-08-02 15:23:36.190587] chlorobium.phaeobacteroides.dsm.266_CP000492_ab: max_subopts=4, avg_subopts=1.12087912088; std_subopts=0.370184534771 chlorobium.phaeobacteroides.dsm.266_CP000492_ab 0.342857142857 [(4, 4, 16, (123, 5))] 2 0.874142665639 0.103323991743 ************************************************** [2017-08-02 15:23:36.259595] scenedesmus.obliquus_DQ396875_n: max_subopts=11, avg_subopts=1.63571428571; std_subopts=1.36561889486 scenedesmus.obliquus_DQ396875_n 0.333333333333 [(1, 59, 3, (47, 10)), (2, 48, 7, (62, 5)), (1, 60, 3, (77, 10)), (2, 47, 7, (47, 5))] 4 0.803470777692 0.156045802626 ************************************************** [2017-08-02 15:23:36.319390] pyrococcus.abyssi.ge5_AJ248288: max_subopts=5, avg_subopts=1.30350877193; std_subopts=0.644127292716 pyrococcus.abyssi.ge5_AJ248288 0.368421052632 [(4, 0, 16, (19, 1))] 2 0.858706788787 0.087924759981 ************************************************** [2017-08-02 15:23:36.386837] euglena.longa_AJ294725_s: max_subopts=6, avg_subopts=1.50308008214; std_subopts=0.832779887251 euglena.longa_AJ294725_s 0.153846153846 [(6, 23, 18, (229, 10))] 3 0.888978717059 0.113922567367 ************************************************** [2017-08-02 15:23:36.438490] marchantia.polymorpha_M.CAU.C.Mar.polym.1: max_subopts=7, avg_subopts=1.20411985019; std_subopts=0.531323935291 marchantia.polymorpha_M.CAU.C.Mar.polym.1 0.341463414634 [(3, 14, 9, (-18, 1))] 2 0.863562399373 0.118631879729 ************************************************** [2017-08-02 15:23:36.473985] lactobacillus.delbrueckii.subsp..bulgaricus.atcc.baa-365_CP000412_cj: max_subopts=8, avg_subopts=1.16597510373; std_subopts=0.591646538007 lactobacillus.delbrueckii.subsp..bulgaricus.atcc.baa-365_CP000412_cj 0.3125 [(6, 7, 21, (236, 5))] 4 0.768513096695 0.158627773623 ************************************************** [2017-08-02 15:23:36.539924] uncultured.alteromonadales.bacterium_AY934393: max_subopts=8, avg_subopts=1.35793357934; std_subopts=0.752070801224 uncultured.alteromonadales.bacterium_AY934393 0.133333333333 [(5, 19, 16, (93, 10))] 3 0.881780775273 0.150398168466 ************************************************** [2017-08-02 15:23:36.573893] methanospirillum.hungatei.jf-1_CP000254_ad: max_subopts=3, avg_subopts=1.09271523179; std_subopts=0.304875613419 methanospirillum.hungatei.jf-1_CP000254_ad 0.4375 [(4, 9, 16, (127, 5))] 3 0.829808774971 0.137522820943 ************************************************** [2017-08-02 15:23:36.632245] drosophila.melanogaster_DROTRF2: max_subopts=8, avg_subopts=1.24773960217; std_subopts=0.662970020928 drosophila.melanogaster_DROTRF2 0.228571428571 [(7, 4, 22, (33, 1))] 2 0.840682945732 0.13748377754 ************************************************** [2017-08-02 15:23:36.696022] pelobacter.carbinolicus.dsm.2380_CP000142_af: max_subopts=6, avg_subopts=1.2785467128; std_subopts=0.644586156818 pelobacter.carbinolicus.dsm.2380_CP000142_af 0.0 [(3, 16, 9, (-159, 5))] 2 0.849592460389 0.130073160151 ************************************************** [2017-08-02 15:23:36.760001] bordetella.parapertussis_BX640433_e: max_subopts=4, avg_subopts=1.41209829868; std_subopts=0.633295713736 bordetella.parapertussis_BX640433_e 0.111111111111 [(7, 8, 21, (211, 10))] 3 0.909841707903 0.125162676024 ************************************************** [2017-08-02 15:23:36.826983] magnetococcus.sp..mc-1_CP000471_ao: max_subopts=8, avg_subopts=1.25985663082; std_subopts=0.654897169047 magnetococcus.sp..mc-1_CP000471_ao 0.0 [(5, 0, 15, (17, 10))] 2 0.798075831352 0.194741307287 ************************************************** [2017-08-02 15:23:36.876035] yarrowia.lipolytica.clib122_CR382129_by: max_subopts=6, avg_subopts=1.36218678815; std_subopts=0.753748870107 yarrowia.lipolytica.clib122_CR382129_by 0.545454545455 [(1, 0, 5, (-179, 10))] 2 0.885014337607 0.0677377810148 ************************************************** [2017-08-02 15:23:36.942848] hyperthermus.butylicus.dsm.5456_CP000493_ad: max_subopts=7, avg_subopts=1.37972508591; std_subopts=0.856630452402 hyperthermus.butylicus.dsm.5456_CP000493_ad 0.266666666667 [(3, 48, 9, (-9, 2))] 4 0.832668707719 0.133543784606 ************************************************** [2017-08-02 15:23:37.003111] homo.sapiens_AC004932_k: max_subopts=8, avg_subopts=1.33333333333; std_subopts=0.834781350643 homo.sapiens_AC004932_k 0.368421052632 [(7, 2, 21, (259, 10))] 3 0.84421928387 0.100033815051 ************************************************** [2017-08-02 15:23:37.068955] uncultured.gamma.proteobacterium_AY934387_a: max_subopts=4, avg_subopts=1.15371024735; std_subopts=0.440099588153 uncultured.gamma.proteobacterium_AY934387_a 0.228571428571 [(5, 13, 15, (-21, 5))] 2 0.824871252527 0.1503357371 ************************************************** [2017-08-02 15:23:37.110242] rickettsia.conorii.str..malish.7_AE008623: max_subopts=4, avg_subopts=1.20212765957; std_subopts=0.50063336978 rickettsia.conorii.str..malish.7_AE008623 0.0952380952381 [(6, 0, 18, (72, 5))] 2 0.817211712968 0.14509640989 ************************************************** [2017-08-02 15:23:37.183170] solibacter.usitatus.ellin6076_CP000473_ay: max_subopts=16, avg_subopts=1.36448598131; std_subopts=1.02120071681 solibacter.usitatus.ellin6076_CP000473_ay 0.428571428571 [(6, 19, 19, (197, 10))] 4 0.90316291976 0.112417787819 ************************************************** [2017-08-02 15:23:37.231498] wigglesworthia.glossinidia.endosymbiont.of.glossina.brevipalpis_BA000021_ad: max_subopts=14, avg_subopts=1.22737306843; std_subopts=0.860695049938 wigglesworthia.glossinidia.endosymbiont.of.glossina.brevipalpis_BA000021_ad 0.1875 [(7, 14, 22, (211, 5))] 6 0.724147866584 0.185507288719 ************************************************** [2017-08-02 15:23:37.267137] sulfolobus.tokodaii.str..7_BA000023_ak: max_subopts=4, avg_subopts=1.14065510597; std_subopts=0.427233285094 sulfolobus.tokodaii.str..7_BA000023_ak 0.341463414634 [(4, 11, 13, (-83, 5))] 2 0.801168633193 0.14612742981 ************************************************** [2017-08-02 15:23:37.332700] uncultured.marinobacter.sp._AY934351: max_subopts=6, avg_subopts=1.31971580817; std_subopts=0.749720083382 uncultured.marinobacter.sp._AY934351 0.473684210526 [(6, 7, 19, (61, 5))] 3 0.853723129365 0.107518972579 ************************************************** [2017-08-02 15:23:37.396886] nephroselmis.olivacea_AF137379_g: max_subopts=6, avg_subopts=1.3368869936; std_subopts=0.75701672222 nephroselmis.olivacea_AF137379_g 0.48275862069 [(4, 8, 17, (224, 5))] 3 0.844223537816 0.0765340074238 **************************************************

for x in stat_accs:
    print x

(0.5384615384615384, 0.83465863237251903, 0.10772401894339267, 0.10772401894339267) (0.22857142857142856, 0.81822699630263351, 0.15525252168871162, 0.15525252168871162) (0.23529411764705882, 0.81093307030350159, 0.14928734299122923, 0.14928734299122923) (0.2727272727272727, 0.80428918510387604, 0.14209353592574656, 0.14209353592574656) (0.2916666666666667, 0.8732875913546192, 0.11183217664142531, 0.11183217664142531) (0.46153846153846156, 0.86320292027278667, 0.12268617275674173, 0.12268617275674173) (0.18181818181818182, 0.84378794615813102, 0.099945414929723309, 0.099945414929723309) (0.1875, 0.76870289329439645, 0.19247994268642413, 0.19247994268642413) (0.32558139534883723, 0.78979415924216112, 0.165393594840964, 0.165393594840964) (0.17647058823529413, 0.80554997467831513, 0.1690767570890182, 0.1690767570890182) (0.16216216216216217, 0.76300147785941541, 0.20142758569635033, 0.20142758569635033) (0.24242424242424243, 0.84056774944787616, 0.15121638928909983, 0.15121638928909983) (0.358974358974359, 0.83799009441112804, 0.13863193474744753, 0.13863193474744753) (0.20512820512820512, 0.84563272842548043, 0.13522223028952163, 0.13522223028952163) (0.28, 0.83463560306728013, 0.14272978548774162, 0.14272978548774162) (0.4444444444444444, 0.84467001322469459, 0.12883433091134094, 0.12883433091134094) (0.21621621621621623, 0.78585637918283258, 0.1771052626828365, 0.1771052626828365) (0.34285714285714286, 0.87414266563946863, 0.10332399174346914, 0.10332399174346914) (0.3333333333333333, 0.80347077769218112, 0.15604580262555093, 0.15604580262555093) (0.3684210526315789, 0.85870678878745654, 0.087924759980970102, 0.087924759980970102) (0.15384615384615385, 0.88897871705932241, 0.11392256736707543, 0.11392256736707543) (0.34146341463414637, 0.86356239937283874, 0.11863187972889236, 0.11863187972889236) (0.3125, 0.76851309669465784, 0.158627773623193, 0.158627773623193) (0.13333333333333333, 0.88178077527281962, 0.15039816846582588, 0.15039816846582588) (0.4375, 0.82980877497074912, 0.13752282094325252, 0.13752282094325252) (0.22857142857142856, 0.84068294573182034, 0.13748377753971389, 0.13748377753971389) (0.0, 0.849592460389054, 0.13007316015114975, 0.13007316015114975) (0.1111111111111111, 0.90984170790275931, 0.12516267602355205, 0.12516267602355205) (0.0, 0.79807583135222115, 0.19474130728687397, 0.19474130728687397) (0.5454545454545454, 0.88501433760721693, 0.067737781014773971, 0.067737781014773971) (0.26666666666666666, 0.83266870771856139, 0.13354378460561847, 0.13354378460561847) (0.3684210526315789, 0.84421928387028811, 0.10003381505137106, 0.10003381505137106) (0.22857142857142856, 0.82487125252746329, 0.15033573709973833, 0.15033573709973833) (0.09523809523809523, 0.81721171296771089, 0.14509640988971242, 0.14509640988971242) (0.42857142857142855, 0.90316291975991358, 0.11241778781855627, 0.11241778781855627) (0.1875, 0.72414786658356523, 0.18550728871863098, 0.18550728871863098) (0.34146341463414637, 0.80116863319308917, 0.14612742980971027, 0.14612742980971027) (0.47368421052631576, 0.85372312936461558, 0.10751897257862499, 0.10751897257862499) (0.4827586206896552, 0.84422353781607029, 0.076534007423825143, 0.076534007423825143)

-

rs = 'AGCUACGUAGCUUAGUUGGUUAGAGCACGGCACUCAUAACGCCGAGGUCACAGGUUCAAUUCCUGUCGUAGCUA'
for st in set([s.replace('>','.').replace('<','.') for s in ss]):
    print ">\n{1}\n{0}".format(st,rs)

> AGCUACGUAGCUUAGUUGGUUAGAGCACGGCACUCAUAACGCCGAGGUCACAGGUUCAAUUCCUGUCGUAGCUA ((...)((.((...)..(((...).(...)).((...).(...).)..).).(...).).(...).(...))). > AGCUACGUAGCUUAGUUGGUUAGAGCACGGCACUCAUAACGCCGAGGUCACAGGUUCAAUUCCUGUCGUAGCUA (((...)(.((...)..(((...).(...)).((...).(...).)..).).(...).).(...).)(...)). > AGCUACGUAGCUUAGUUGGUUAGAGCACGGCACUCAUAACGCCGAGGUCACAGGUUCAAUUCCUGUCGUAGCUA ((...)((.((...)..(((...).(...)).((...).(...).)..).).(...).).(...).)(...)). > AGCUACGUAGCUUAGUUGGUUAGAGCACGGCACUCAUAACGCCGAGGUCACAGGUUCAAUUCCUGUCGUAGCUA (.(.((...)(...).)(.(...)((.((...)(...).(...).)(...).).(...)..)).).(...).). > AGCUACGUAGCUUAGUUGGUUAGAGCACGGCACUCAUAACGCCGAGGUCACAGGUUCAAUUCCUGUCGUAGCUA (.(.((...)(...).)(.(...)((..(...).(.(...)(...)(...).))(...)..)).).(...).). > AGCUACGUAGCUUAGUUGGUUAGAGCACGGCACUCAUAACGCCGAGGUCACAGGUUCAAUUCCUGUCGUAGCUA (.(.((...)(...).)(.(...)((...)(.((...).(...).)(...).).(...)..)).).(...).).

print minsig

[(6, 25, 18, (133, 10))]

`list(combinations([1,2],2))

[(1, 2)]

list(combinations(['<(><(>(...)>(....)>(...)(...))><(...)><((...)((...)(...)(...)>.<(...)>)>)>)>'],2))

[]

names

['uncultured.alpha.proteobacterium_AY933788', 'magnetococcus.sp..mc-1_CP000471_ao', 'uncultured.thiotrichales.bacterium_AY934254_a', 'carboxydothermus.hydrogenoformans.z-2901_CP000141_af', 'nocardia.farcinica.ifm.10152_AP006618_az', 'hyperthermus.butylicus.dsm.5456_CP000493_ad', 'methanospirillum.hungatei.jf-1_CP000254_ad', 'synechococcus.sp..wh.8102_BX569691', 'wigglesworthia.glossinidia.endosymbiont.of.glossina.brevipalpis_BA000021_ad', 'yersinia.pestis.kim_AE009952_ae', 'uncultured.alpha.proteobacterium_DQ093144_a', 'uncultured.alteromonadales.bacterium_AY934393', 'pyrococcus.abyssi.ge5_AJ248288', 'drosophila.melanogaster_DROTRF2', 'arabidopsis.thaliana_AC006340', 'lactobacillus.casei.atcc.334_CP000423_j', 'uncultured.gamma.proteobacterium_AY934387_a', 'sulfolobus.tokodaii.str..7_BA000023_ak', 'euglena.longa_AJ294725_s', 'methanosarcina.barkeri.str..fusaro_CP000099_bi', 'chlorobium.phaeobacteroides.dsm.266_CP000492_ab', 'rickettsia.conorii.str..malish.7_AE008623', 'solibacter.usitatus.ellin6076_CP000473_ay', 'homo.sapiens_X04779', 'bordetella.parapertussis_BX640433_e', 'uncultured.gamma.proteobacterium_AY934184_a', 'tropheryma.whipplei.str..twist_AE014184_h', 'nephroselmis.olivacea_AF137379_g', 'chlamydophila.pneumoniae.tw-183_AE017159_j', 'haloferax.volcanii_HALTRK1', 'thermoplasma.volcanium.gss1_BA000011_b', 'pelobacter.carbinolicus.dsm.2380_CP000142_af', 'yarrowia.lipolytica.clib122_CR382129_by', 'homo.sapiens_AC004932_k', 'uncultured.marinobacter.sp._AY934351', 'marchantia.polymorpha_M.CAU.C.Mar.polym.1', 'scenedesmus.obliquus_DQ396875_n', 'nitrosomonas.europaea.atcc.19718_BX321863_b', 'lactobacillus.delbrueckii.subsp..bulgaricus.atcc.baa-365_CP000412_cj']

t = timestamp('')
# 3 minutes
classic_slices,classic_accuracies = classic_slices_and_accuracy(names)
print timestamp(t)

[2017-09-14 17:37:21.680657] [2017-09-14 17:33:23.312065]

from numpy import average as avg
accuracies = []
for k in classic_accuracies:
    accuracies.append(avg(classic_accuracies[k]))
print avg(accuracies)
print std(accuracies)
print min(accuracies), len([x for x in accuracies if x==min(accuracies)])
print max(accuracies), len([x for x in accuracies if x==max(accuracies)])
print med(accuracies)
print accuracies

0.568943997152 0.268530539034 0.0 1 1.0 2 0.553191489362 [0.96399535423925664, 0.83282828282828292, 0.21276595744680851, 0.69767441860465118, 0.24390243902439024, 0.31111111111111112, 0.58536585365853655, 0.55319148936170215, 1.0, 0.0, 0.77272727272727271, 0.4642857142857143, 0.76190476190476186, 0.80000000000000004, 0.71111111111111114, 0.25, 0.32558139534883723, 0.30434782608695654, 0.85365853658536583, 0.22727272727272727, 0.69565217391304346, 1.0, 0.20782726045883937, 0.30434782608695654, 0.58434343434343439, 0.2978723404255319, 0.93841463414634152, 0.52753623188405796, 0.34999999999999998, 0.46511627906976744, 0.8290493676050591, 0.25641025641025639, 0.72727272727272729, 0.51618871415356149, 0.54545454545454541, 0.5, 0.97674418604651159, 0.91304347826086951, 0.68181818181818177]

%md
# Classic slice accuracy
The average accuracy for the classic parameters (considering suboptimal structures) is 0.568943997152 (median 0.553191489362) with std 0.268530539034.

Classic slice accuracy

The average accuracy for the classic parameters (considering suboptimal structures) is 0.568943997152 (median 0.553191489362) with std 0.268530539034.

common_pt = (489/40, 51/320, -231/80)

t=timestamp('')
# 3 minutes
mavg_slices,tweaked_accuracies = pt_slices_and_accuracy(names,common_pt)
print timestamp(t)

[2017-09-14 17:43:02.699808] [2017-09-14 17:38:03.732423]

accuracies = []
for k in tweaked_accuracies:
    accuracies.append(avg(tweaked_accuracies[k]))
print avg(accuracies)
print std(accuracies)
print min(accuracies), len([x for x in accuracies if x==min(accuracies)])
print max(accuracies), len([x for x in accuracies if x==max(accuracies)])
print med(accuracies)
print accuracies

0.742401994833 0.252975709202 0.0 1 1.0 6 0.772727272727 [0.96399535423925664, 0.71111111111111114, 0.21276595744680851, 0.69767441860465118, 1.0, 0.97674418604651159, 0.55813953488372092, 0.78260869565217395, 1.0, 0.0, 0.77272727272727271, 0.5, 0.48888888888888887, 0.80000000000000004, 0.71111111111111114, 0.68085106382978722, 0.97560975609756095, 0.97674418604651159, 0.85365853658536583, 0.68292682926829273, 0.69565217391304346, 1.0, 0.20782726045883937, 0.97674418604651159, 1.0, 0.78260869565217395, 0.93841463414634152, 0.96564482029598309, 0.75, 1.0, 0.69767441860465118, 0.24390243902439024, 0.72727272727272729, 0.50531914893617014, 0.54545454545454541, 1.0, 0.97674418604651159, 0.91304347826086951, 0.68181818181818177]

%md
# Tweaking the parameters
The average accuracy after tweaking the parameters (considering suboptimal structures) is 0.742401994833 (median 0.772727272727) with std 0.252975709202.

Tweaking the parameters

The average accuracy after tweaking the parameters (considering suboptimal structures) is 0.742401994833 (median 0.772727272727) with std 0.252975709202.

t = timestamp('')
# 5 min
classic_matrices,classic_bounded_matrices=similarity_matrices(names,classic_pt)
print timestamp(t)

[2017-07-31 04:01:08.313040] 0 Avoiding self [2017-07-31 04:01:12.849542] 1 Avoiding self [2017-07-31 04:01:16.368548] 2 Avoiding self [2017-07-31 04:01:40.872290] 3 Avoiding self [2017-07-31 04:01:46.979064] 4 Avoiding self [2017-07-31 04:01:52.453668] 5 Avoiding self [2017-07-31 04:01:57.591340] 6 Avoiding self [2017-07-31 04:02:02.632839] 7 Avoiding self [2017-07-31 04:02:08.367036] 8 Avoiding self Avoiding self [2017-07-31 04:02:23.054167] 9 Avoiding self [2017-07-31 04:02:29.800010] 10 Avoiding self [2017-07-31 04:02:43.141657] 11 Avoiding self [2017-07-31 04:02:50.896981] 12 Avoiding self [2017-07-31 04:02:54.528857] 13 Avoiding self [2017-07-31 04:03:09.158544] 14 Avoiding self [2017-07-31 04:03:22.011483] 15 Avoiding self [2017-07-31 04:03:36.224881] 16 Avoiding self Avoiding self [2017-07-31 04:03:46.146898] 17 Avoiding self [2017-07-31 04:03:49.379607] 18 Avoiding self [2017-07-31 04:03:52.300318] 19 Avoiding self [2017-07-31 04:04:01.229156] 20 Avoiding self [2017-07-31 04:04:04.777856] 21 Avoiding self [2017-07-31 04:04:08.601774] 22 Avoiding self [2017-07-31 04:04:16.899862] 23 Avoiding self [2017-07-31 04:04:20.280625] 24 Avoiding self [2017-07-31 04:04:30.101112] 25 Avoiding self [2017-07-31 04:04:34.520678] 26 Avoiding self [2017-07-31 04:04:39.147908] 27 Avoiding self [2017-07-31 04:04:49.454179] 28 Avoiding self [2017-07-31 04:05:02.107538] 29 Avoiding self Avoiding self [2017-07-31 04:05:06.510930] 30 Avoiding self [2017-07-31 04:05:10.424936] 31 Avoiding self [2017-07-31 04:05:21.042569] 32 Avoiding self [2017-07-31 04:05:25.112133] 33 Avoiding self [2017-07-31 04:05:29.397879] 34 Avoiding self [2017-07-31 04:05:41.525826] 35 Avoiding self Avoiding self Avoiding self [2017-07-31 04:05:54.155492] 36 Avoiding self [2017-07-31 04:06:04.290441] 37 Avoiding self [2017-07-31 04:06:16.913338] 38 Avoiding self [2017-07-31 04:06:20.357623] [2017-07-31 04:01:08.312665]

def analyze_matrices(classic_matrices):
    mins = []
    maxs = []
    avgs = []
    othr = []
    for name in classic_matrices:
        slice_mins = []
        slice_avgs = []
        slice_maxs = []
        slice_othr = []
        for sl in classic_matrices[name]:
            try:
                min_nbr = min([min([min(row) for row in classic_matrices[name][sl][nbr]]) for nbr in classic_matrices[name][sl]])
            except ValueError:
                continue
            max_nbr = max([max([max(row) for row in classic_matrices[name][sl][nbr]]) for nbr in classic_matrices[name][sl]])
            avg_nbr = avg([avg([avg(row) for row in classic_matrices[name][sl][nbr]]) for nbr in classic_matrices[name][sl]])
            oth_nbr = min([min([max(row) for row in classic_matrices[name][sl][nbr]]) for nbr in classic_matrices[name][sl]])

            slice_mins.append(min_nbr)
            slice_maxs.append(max_nbr)
            slice_avgs.append(avg_nbr)
            slice_othr.append(oth_nbr)

        try:
            mins.append(min(slice_mins))
        except ValueError:
            continue
        maxs.append(max(slice_maxs))
        avgs.append(avg(slice_avgs))
        #if len(slice_othr)>1:
        #    continue
        othr.extend(slice_othr)

    return mins, avgs, maxs, othr

cmins,cavgs,cmaxs,co=analyze_matrices(classic_matrices)

avg(co)
std(co)

0.29678980258415633 0.20060049060704496

sorted(co)

[0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0425531914893617, 0.05405405405405406, 0.15, 0.20408163265306123, 0.2127659574468085, 0.2222222222222222, 0.24489795918367346, 0.25, 0.2553191489361702, 0.26666666666666666, 0.2857142857142857, 0.2857142857142857, 0.2857142857142857, 0.3, 0.3181818181818182, 0.32558139534883723, 0.35, 0.35555555555555557, 0.4444444444444444, 0.4583333333333333, 0.46511627906976744, 0.46511627906976744, 0.4782608695652174, 0.4878048780487805, 0.4878048780487805, 0.5, 0.5, 0.5116279069767442, 0.5116279069767442, 0.52, 0.5238095238095238, 0.55, 0.5581395348837209, 0.5789473684210527, 0.6086956521739131]

%md
# Classic Robustness
Using the structures in the classic slice as reference, it is expected that any structure from a neighboring region will have accuracy 0.3 with respect with respect to at least one structure from the classic slice.

Classic Robustness

Using the structures in the classic slice as reference, it is expected that any structure from a neighboring region will have accuracy 0.3 with respect with respect to at least one structure from the classic slice.

0.13028404049096914

|max_common = 0
max_names = []
for n in classic_matrices:
    if len(classic_matrices[n])>=2:
        #if max_common == len(classic_matrices[n]):
        max_names.append(n)
        #else:
        #    max_common = len(classic_matrices[n])
        #    max_names=[n]
print max_names

['uncultured.alpha.proteobacterium_DQ093144_a', 'nocardia.farcinica.ifm.10152_AP006618_az', 'yarrowia.lipolytica.clib122_CR382129_by', 'wigglesworthia.glossinidia.endosymbiont.of.glossina.brevipalpis_BA000021_ad']

%md
The classic parameters point is in 3 regions for wigglesworthia.glossinidia.endosymbiont.of.glossina.brevipalpis_BA000021_ad

It belongs to 2 regions for uncultured.alpha.proteobacterium_DQ093144_a, nocardia.farcinica.ifm.10152_AP006618_az, andyarrowia.lipolytica.clib122_CR382129_by

The classic parameters point is in 3 regions for wigglesworthia.glossinidia.endosymbiont.of.glossina.brevipalpis_BA000021_ad

It belongs to 2 regions for uncultured.alpha.proteobacterium_DQ093144_a, nocardia.farcinica.ifm.10152_AP006618_az, andyarrowia.lipolytica.clib122_CR382129_by

t = timestamp('')
tweaked_matrices,tweaked_bounded_matrices=similarity_matrices(names,common_pt)
print timestamp(t)

[2017-07-31 04:24:24.685636] 0 Avoiding self [2017-07-31 04:24:35.130766] 1 Avoiding self [2017-07-31 04:24:39.796842] 2 Avoiding self [2017-07-31 04:24:44.898496] 3 Avoiding self [2017-07-31 04:24:50.634261] 4 Avoiding self [2017-07-31 04:24:56.706563] 5 Avoiding self [2017-07-31 04:25:01.494707] 6 Avoiding self [2017-07-31 04:25:06.213420] 7 Avoiding self [2017-07-31 04:25:11.370366] 8 Avoiding self [2017-07-31 04:25:17.002026] 9 Avoiding self [2017-07-31 04:25:23.160244] 10 Avoiding self [2017-07-31 04:25:26.987425] 11 Avoiding self [2017-07-31 04:25:32.624599] 12 Avoiding self [2017-07-31 04:25:38.566758] 13 Avoiding self [2017-07-31 04:25:41.321708] 14 Avoiding self [2017-07-31 04:25:48.050045] 15 Avoiding self [2017-07-31 04:25:58.385046] 16 Avoiding self [2017-07-31 04:26:05.933164] 17 Avoiding self [2017-07-31 04:26:09.056931] 18 Avoiding self [2017-07-31 04:26:15.752849] 19 Avoiding self [2017-07-31 04:26:19.977103] 20 Avoiding self [2017-07-31 04:26:23.703478] 21 Avoiding self [2017-07-31 04:26:31.498553] 22 Avoiding self [2017-07-31 04:26:35.066069] 23 Avoiding self [2017-07-31 04:26:43.422917] 24 Avoiding self [2017-07-31 04:26:47.461248] 25 Avoiding self [2017-07-31 04:26:51.592168] 26 Avoiding self [2017-07-31 04:27:00.928606] 27 Avoiding self [2017-07-31 04:27:04.605150] 28 Avoiding self [2017-07-31 04:27:08.237030] 29 Avoiding self [2017-07-31 04:27:11.161391] 30 Avoiding self [2017-07-31 04:27:20.850958] 31 Avoiding self [2017-07-31 04:27:24.968905] 32 Avoiding self [2017-07-31 04:27:29.315101] 33 Avoiding self [2017-07-31 04:27:39.673430] 34 Avoiding self [2017-07-31 04:27:43.592278] 35 Avoiding self [2017-07-31 04:27:46.572899] 36 Avoiding self [2017-07-31 04:27:50.339650] 37 Avoiding self [2017-07-31 04:27:54.303687] 38 Avoiding self [2017-07-31 04:28:04.805764] [2017-07-31 04:24:24.683613]

twmins,twavgs,twmaxs,two=analyze_matrices(tweaked_matrices)

sorted(two)

[0.0, 0.0, 0.0, 0.20408163265306123, 0.2127659574468085, 0.2127659574468085, 0.21818181818181817, 0.23809523809523808, 0.24489795918367346, 0.25, 0.25, 0.2608695652173913, 0.2692307692307692, 0.2727272727272727, 0.2857142857142857, 0.2916666666666667, 0.2926829268292683, 0.2978723404255319, 0.3, 0.3181818181818182, 0.32558139534883723, 0.35555555555555557, 0.35555555555555557, 0.3684210526315789, 0.4444444444444444, 0.4583333333333333, 0.46511627906976744, 0.46511627906976744, 0.46808510638297873, 0.4878048780487805, 0.4897959183673469, 0.5, 0.5, 0.5, 0.5116279069767442, 0.52, 0.5581395348837209, 0.5714285714285714, 0.6086956521739131]

avg(two)
std(two)

0.34290860695567454 0.15236304375062992

# 3.5 minutes
data_dict = {}
n = len(sobj_names)
for i,sobj_name in enumerate(sobj_names):
    name = names[i]
    print timestamp("[{0: 3}/{1}; {2: 5.2f}] Started {3}".format(i+1,n,float(100*i/float(n)),name))
    seq_dict = {}

    P = RNAPolytope.construct_from_file(sobj_name)
    classic_slices = get_classic_slices(P)
    seq_dict['classic_pt_intervals'] = []
    seq_dict['classic_rg_intervals'] = []
    pt_sensitivity = Infinity
    rg_sensitivity = Infinity

    for sl in classic_slices:

        intervals = get_var_intervals(sl,'classic')
        seq_dict['classic_pt_intervals'].append(intervals)
        if min(intervals)< pt_sensitivity:
            pt_sensitivity = min(intervals)


        intervals = get_var_intervals(sl)
        seq_dict['classic_rg_intervals'].append(intervals)
        if min(intervals)< rg_sensitivity:
            rg_sensitivity = min(intervals)

    seq_dict['classic_pt_sensitivity']= pt_sensitivity
    seq_dict['classic_rg_sensitivity']= rg_sensitivity

    max_avg_slices = get_max_average_slice(P,name)
    seq_dict['mavg_intervals'] = []
    mavg_sensitivity = Infinity

    for sl in max_avg_slices:

        intervals = get_var_intervals(sl)
        seq_dict['mavg_intervals'].append(intervals)
        if min(intervals)< mavg_sensitivity:
            mavg_sensitivity = min(intervals)

    seq_dict['mavg_sensitivity'] = mavg_sensitivity
    data_dict[name] = seq_dict

[2017-07-30 19:38:30.807540] [ 1/39; 0.00] Started nitrosomonas.europaea.atcc.19718_BX321863_b [2017-07-30 19:38:33.653337] [ 2/39; 2.56] Started lactobacillus.casei.atcc.334_CP000423_j [2017-07-30 19:38:36.000858] [ 3/39; 5.13] Started methanosarcina.barkeri.str..fusaro_CP000099_bi [2017-07-30 19:38:38.628835] [ 4/39; 7.69] Started uncultured.alpha.proteobacterium_AY933788 [2017-07-30 19:38:41.265167] [ 5/39; 10.26] Started tropheryma.whipplei.str..twist_AE014184_h [2017-07-30 19:38:58.575148] [ 6/39; 12.82] Started uncultured.thiotrichales.bacterium_AY934254_a [2017-07-30 19:39:02.886653] [ 7/39; 15.38] Started synechococcus.sp..wh.8102_BX569691 [2017-07-30 19:39:06.658890] [ 8/39; 17.95] Started yersinia.pestis.kim_AE009952_ae [2017-07-30 19:39:11.036377] [ 9/39; 20.51] Started uncultured.alpha.proteobacterium_DQ093144_a [2017-07-30 19:39:15.557219] [ 10/39; 23.08] Started carboxydothermus.hydrogenoformans.z-2901_CP000141_af [2017-07-30 19:39:20.094671] [ 11/39; 25.64] Started haloferax.volcanii_HALTRK1 [2017-07-30 19:39:23.995749] [ 12/39; 28.21] Started homo.sapiens_X04779 [2017-07-30 19:39:28.562320] [ 13/39; 30.77] Started arabidopsis.thaliana_AC006340 [2017-07-30 19:39:32.356698] [ 14/39; 33.33] Started chlamydophila.pneumoniae.tw-183_AE017159_j [2017-07-30 19:39:36.424664] [ 15/39; 35.90] Started thermoplasma.volcanium.gss1_BA000011_b [2017-07-30 19:39:40.849930] [ 16/39; 38.46] Started uncultured.gamma.proteobacterium_AY934184_a [2017-07-30 19:39:44.991796] [ 17/39; 41.03] Started nocardia.farcinica.ifm.10152_AP006618_az [2017-07-30 19:39:49.789360] [ 18/39; 43.59] Started chlorobium.phaeobacteroides.dsm.266_CP000492_ab [2017-07-30 19:39:53.617117] [ 19/39; 46.15] Started scenedesmus.obliquus_DQ396875_n [2017-07-30 19:39:56.698762] [ 20/39; 48.72] Started pyrococcus.abyssi.ge5_AJ248288 [2017-07-30 19:40:00.845546] [ 21/39; 51.28] Started euglena.longa_AJ294725_s [2017-07-30 19:40:04.790671] [ 22/39; 53.85] Started marchantia.polymorpha_M.CAU.C.Mar.polym.1 [2017-07-30 19:40:09.033923] [ 23/39; 56.41] Started lactobacillus.delbrueckii.subsp..bulgaricus.atcc.baa-365_CP000412_cj [2017-07-30 19:40:13.184600] [ 24/39; 58.97] Started uncultured.alteromonadales.bacterium_AY934393 [2017-07-30 19:40:16.832205] [ 25/39; 61.54] Started methanospirillum.hungatei.jf-1_CP000254_ad [2017-07-30 19:40:20.937581] [ 26/39; 64.10] Started drosophila.melanogaster_DROTRF2 [2017-07-30 19:40:24.859138] [ 27/39; 66.67] Started pelobacter.carbinolicus.dsm.2380_CP000142_af [2017-07-30 19:40:29.362620] [ 28/39; 69.23] Started bordetella.parapertussis_BX640433_e [2017-07-30 19:40:33.623474] [ 29/39; 71.79] Started magnetococcus.sp..mc-1_CP000471_ao [2017-07-30 19:40:38.173285] [ 30/39; 74.36] Started yarrowia.lipolytica.clib122_CR382129_by [2017-07-30 19:40:42.003734] [ 31/39; 76.92] Started hyperthermus.butylicus.dsm.5456_CP000493_ad [2017-07-30 19:40:46.331449] [ 32/39; 79.49] Started homo.sapiens_AC004932_k [2017-07-30 19:40:50.494343] [ 33/39; 82.05] Started uncultured.gamma.proteobacterium_AY934387_a [2017-07-30 19:40:55.035128] [ 34/39; 84.62] Started rickettsia.conorii.str..malish.7_AE008623 [2017-07-30 19:40:58.489881] [ 35/39; 87.18] Started solibacter.usitatus.ellin6076_CP000473_ay [2017-07-30 19:41:03.415520] [ 36/39; 89.74] Started wigglesworthia.glossinidia.endosymbiont.of.glossina.brevipalpis_BA000021_ad [2017-07-30 19:41:07.298645] [ 37/39; 92.31] Started sulfolobus.tokodaii.str..7_BA000023_ak [2017-07-30 19:41:11.470772] [ 38/39; 94.87] Started uncultured.marinobacter.sp._AY934351 [2017-07-30 19:41:16.277704] [ 39/39; 97.44] Started nephroselmis.olivacea_AF137379_g

from numpy import average as avg, float32 as npfloat
ab_ratios=[]
ac_ratios=[]
bc_ratios=[]
for n in data_dict:
    for a,b,c in data_dict[n]['classic_pt_intervals']:
        if b!=0:
            ab_ratios.append(npfloat(a/b))
        if c!=0:
            ac_ratios.append(npfloat(a/c))
            bc_ratios.append(npfloat(b/c))
print "Classic pt intervals"
print "AB Ratios"
avg(ab_ratios)
std(ab_ratios)
print "-"*10

print "AC Ratios"
avg(ac_ratios)
std(ac_ratios)
print "-"*10

print "BC Ratios"
avg(bc_ratios)
std(bc_ratios)
print "-"*10
ab_ratios=[]
ac_ratios=[]
bc_ratios=[]

for n in data_dict:
    for a,b,c in data_dict[n]['classic_rg_intervals']:
        if b!=0:
            ab_ratios.append(npfloat(a/b))
        if c!=0:
            ac_ratios.append(npfloat(a/c))
            bc_ratios.append(npfloat(b/c))
print("*"*40)
print "Classic rg intervals"
print "AB Ratios"
avg(ab_ratios)
std(ab_ratios)
print "-"*10

print "AC Ratios"
avg(ac_ratios)
std(ac_ratios)
print "-"*10

print "BC Ratios"
avg(bc_ratios)
std(bc_ratios)
print "-"*10

Classic pt intervals AB Ratios 29.636765 35.3841594242465 ---------- AC Ratios 6.0702877 6.745306689622802 ---------- BC Ratios 0.31000602 0.39570656095909434 ---------- **************************************** Classic rg intervals AB Ratios 67.575783 60.948507465666744 ---------- AC Ratios 4.6494031 2.179624810918347 ---------- BC Ratios 0.19813171 0.5123017501107 ----------

ab_ratios=[]
ac_ratios=[]
bc_ratios=[]

for n in data_dict:
    for a,b,c in data_dict[n]['mavg_intervals']:
        if b!=0:
            ab_ratios.append(npfloat(a/b))
        if c!=0:
            ac_ratios.append(npfloat(a/c))
            bc_ratios.append(npfloat(b/c))
print "Max average"
print "AB Ratios"
avg(ab_ratios)
std(ab_ratios)
print "-"*10

print "AC Ratios"
avg(ac_ratios)
std(ac_ratios)
print "-"*10

print "BC Ratios"
avg(bc_ratios)
std(bc_ratios)
print "-"*10

Max average AB Ratios 47.278629 78.92056424082578 ---------- AC Ratios 5.547956 4.9745379067044295 ---------- BC Ratios 1.4274589 4.95978512627375 ----------

for n in data_dict:
    print data_dict[n]
    break

{'classic_pt_intervals': [[16/5, 1/30, 7/20]], 'mavg_intervals': [[69551/7040, 1/10, 2650643/1182720]], 'mavg_sensitivity': 1/10, 'classic_rg_sensitivity': 5383007891/2732730000, 'classic_pt_sensitivity': 1/30, 'classic_rg_intervals': [[36060573539/327927600, 5383007891/2732730000, 9896110933/609008400]]}

for n in data_dict:
    print data_dict[n]
    break

{'classic_pt_intervals': [[33/5, 1/30, 3/4]], 'mavg_intervals': [[2024749/295680, 4183/29568, 2024749/1182720]], 'mavg_sensitivity': 4183/29568, 'classic_rg_sensitivity': 44835799471/16396380000, 'classic_pt_sensitivity': 1/30, 'classic_rg_intervals': [[16137699043/65585520, 44835799471/16396380000, 5199128197/170522352]]}

from numpy import average as avg, std, median as med
cl_rg_avg = 0
cl_pt_avg = 0

cl_rg_s = [min([data_dict[n]['classic_rg_intervals'][0][0],data_dict[n]['classic_rg_intervals'][0][2]]) for n in data_dict]
cl_pt_s = [min([data_dict[n]['classic_pt_intervals'][0][0],data_dict[n]['classic_pt_intervals'][0][2]]) for n in data_dict]
cl_pt_s = [data_dict[n]['classic_pt_sensitivity'] for n in data_dict]
s_pt = sorted(cl_pt_s)
print s_pt
s_pt = s_pt
print s_pt[4:]
print sum(s_pt)
print std(s_pt)/len(s_pt)
print avg(s_pt)
print med(s_pt)
#print sorted(cl_rg_s)
#print avg(cl_rg_s)
#print std(cl_rg_s)
#print med(cl_rg_s)

[0, 0, 0, 0, 1/70, 1/50, 1/40, 1/40, 1/40, 1/40, 1/30, 1/25, 1/25, 1/20, 1/20, 1/20, 1/20, 1/15, 4/55, 3/40, 7/80, 1/10, 1/10, 1/10, 11/90, 11/90, 1/8, 3/20, 23/140, 1/5, 37/180, 17/80, 13/60, 19/85, 21/85, 1/4, 14/45, 17/50, 29/85] [1/70, 1/50, 1/40, 1/40, 1/40, 1/40, 1/30, 1/25, 1/25, 1/20, 1/20, 1/20, 1/20, 1/15, 4/55, 3/40, 7/80, 1/10, 1/10, 1/10, 11/90, 11/90, 1/8, 3/20, 23/140, 1/5, 37/180, 17/80, 13/60, 19/85, 21/85, 1/4, 14/45, 17/50, 29/85] 387943/3534300 0.109765158589 0.075 0.0976890697461

float(107/222)

0.481981981981982

float(23/168)
float(1/40)

0.13690476190476192 0.025

G

Error in lines 1-1 Traceback (most recent call last): File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/smc_sagews/sage_server.py", line 995, in execute exec compile(block+'\n', '', 'single') in namespace, locals File "", line 1, in <module> NameError: name 'G' is not defined

float(1/(7))

0.14285714285714285

tst = timestamp("")
# 5 minutes
slices,slice2name = get_max_average_slices()
print timestamp(tst)

[2017-08-21 01:04:03.175571] [2017-08-21 01:01:22.636278]

tst = timestamp("")
G,l2s = max_avg_intersection_graph(list(chain(*slices.values())))
print timestamp(tst)

[2017-08-21 01:06:26.208813] [2017-08-21 01:05:00.258245]

G.num_verts()

39

cl=greedy_clique_partition(G)
len(cl)
for c in cl:
    print len(c)

14 15 4 4 3 2 2 2 1 1 1 1 1 1 1

for x in [l2s[x] for x in cl[0]]:
    print "{0}: {1}".format(x.original_vertex.vector(),slice2name[x])

(1, 8, 4, -121/5): euglena.longa_AJ294725_s (1, 8, 4, -303/10): homo.sapiens_X04779 (1, 8, 4, -323/10): uncultured.marinobacter.sp._AY934351 (1, 8, 4, -178/5): nitrosomonas.europaea.atcc.19718_BX321863_b (1, 8, 4, -271/10): uncultured.alteromonadales.bacterium_AY934393 (1, 10, 4, -159/5): lactobacillus.casei.atcc.334_CP000423_j (1, 3, 4, -193/5): haloferax.volcanii_HALTRK1 (1, 8, 4, -309/10): rickettsia.conorii.str..malish.7_AE008623 (1, 8, 4, -179/5): thermoplasma.volcanium.gss1_BA000011_b (1, 5, 4, -291/10): yarrowia.lipolytica.clib122_CR382129_by (1, 9, 4, -157/5): bordetella.parapertussis_BX640433_e (1, 8, 4, -369/10): solibacter.usitatus.ellin6076_CP000473_ay (1, 8, 4, -144/5): uncultured.alpha.proteobacterium_AY933788 (1, 8, 4, -229/10): marchantia.polymorpha_M.CAU.C.Mar.polym.1 (1, 8, 4, -207/5): pelobacter.carbinolicus.dsm.2380_CP000142_af

slices = [ l2s[v] for v in cl[0]]
print(len(slices))

15

slices.sort(key=lambda x: float(x.volume()))

for sl in slices:
    print float(sl.volume())

1.60216666667 2.79823809524 6.32125 18.5400521881 29.1843079167 33.9231887755 34.1336042529 36.4903839323 39.8010432099 50.0135453492 66.233499033 67.5312682427 82.0087287671 100.560111567 140.699411343

def intersect_slices(slices):
    inter = slices[0]
    for sl in slices[1:]:
        inter = inter.intersection(sl)
    return inter
common=intersect_slices(slices)
print common.vertices()
print common.volume()
print get_var_intervals(common)
print [float(x) for x in get_var_intervals(common)]
print float(common.volume())
print [float(x) for x in slice_center(common)]

(A vertex at (129/10, 3/20, -3), A vertex at (123/10, 3/20, -3), A vertex at (12, 3/16, -57/20), A vertex at (117/10, 3/20, -27/10)) 9/8000 [3/20, 3/160, 3/80] [0.15, 0.01875, 0.0375] 0.001125 [12.225, 0.159375, -2.8875]

inter = Polyhedron(((129/10, 3/20, -3),(123/10, 3/20, -3),(12, 3/16, -57/20),(117/10, 3/20, -27/10)))
inter.volume()
slice_center(inter)
[float(x) for x in slice_center(inter)]

9/8000 (489/40, 51/320, -231/80) [12.225, 0.159375, -2.8875]

def intersect_slices(slices):
    inter = slices[0]
    for sl in slices[1:]:
        inter = inter.intersection(sl)
    return inter

for sl in slices:
    print "{0}".format(slice2name[sl])
    print "Volume: {0: 8.2f}; center: ({1: 8.2f},{2: 8.2f},{3: 8.2f})".format(float(sl.volume()),float(sl.center()[0]),float(sl.center()[1]),float(sl.center()[2]))

print "*"*40
from itertools import combinations
mi = 2
max_vol = [0 for i in range(mi+2)]
avoid_list = [[] for i in range(mi+2)]
for i in range(1,mi+2):
    print timestamp(str(i))
    for avoid in combinations(range(len(slices)),i):

        remaining = [j for j in range(len(slices)) if j not in avoid]

        inter = intersect_slices([slices[j] for j in remaining])

        if inter.volume()>max_vol[i]:
            max_vol[i]=inter.volume()
            avoid_list[i] = avoid
            if max_vol[i]>max_vol[0]:
                max_vol[0]=max_vol[i]
                avoid_list[0] = avoid_list[i]

marchantia.polymorpha_M.CAU.C.Mar.polym.1 Volume: 1.60; center: ( 11.83, 0.18, -3.24) yarrowia.lipolytica.clib122_CR382129_by Volume: 2.80; center: ( 14.48, 0.24, -3.24) lactobacillus.casei.atcc.334_CP000423_j Volume: 6.32; center: ( 9.85, 0.01, -2.34) uncultured.marinobacter.sp._AY934351 Volume: 18.54; center: ( 10.30, -0.07, -2.71) uncultured.alpha.proteobacterium_AY933788 Volume: 29.18; center: ( 11.74, -0.02, -3.05) bordetella.parapertussis_BX640433_e Volume: 33.92; center: ( 13.39, -0.09, -3.05) uncultured.alteromonadales.bacterium_AY934393 Volume: 34.13; center: ( 9.76, -0.04, -2.41) nitrosomonas.europaea.atcc.19718_BX321863_b Volume: 36.49; center: ( 12.05, 0.13, -3.27) homo.sapiens_X04779 Volume: 39.80; center: ( 10.17, 0.31, -2.89) thermoplasma.volcanium.gss1_BA000011_b Volume: 50.01; center: ( 13.04, -0.17, -4.46) solibacter.usitatus.ellin6076_CP000473_ay Volume: 66.23; center: ( 14.53, -0.12, -4.81) pelobacter.carbinolicus.dsm.2380_CP000142_af Volume: 67.53; center: ( 10.30, 0.07, -2.95) euglena.longa_AJ294725_s Volume: 82.01; center: ( 8.44, 0.12, -1.80) rickettsia.conorii.str..malish.7_AE008623 Volume: 100.56; center: ( 16.13, -0.18, -4.94) haloferax.volcanii_HALTRK1 Volume: 140.70; center: ( 15.60, 0.83, -6.71) **************************************** [2017-07-29 11:05:58.790837] 1 [2017-07-29 11:06:04.351930] 2 [2017-07-29 11:06:42.619366] 3

max_vol
[float(x) for x in max_vol]

[83717/36000, 1267/14400, 587/1350, 83717/36000] [2.3254722222222224, 0.08798611111111111, 0.4348148148148148, 2.3254722222222224]

avoid_list

[(0, 1, 14), (2,), (1, 2), (0, 1, 14)]

common1 = intersect_slices([sl for i,sl in enumerate(slices) if i not in (0,1,14)])
common1.volume()
float(common1.volume())
common1
get_var_intervals(common1)
[float(x) for x in get_var_intervals(common1)]
slice_center(common1)

83717/36000 2.3254722222222224 A 3-dimensional polyhedron in QQ^3 defined as the convex hull of 13 vertices [1429/780, 19/780, 1429/3120] [1.832051282051282, 0.02435897435897436, 0.4580128205128205] (1381/156, 59/780, -448/195)

print [float(x) for x in get_var_intervals(sl)]
print float(sl.volume())

[0.15, 0.01875, 0.0375] 0.001125

float(sl.volume())

0.34897222222222224

sobj_name='sobjs/uncultured.thiotrichales.bacterium_AY934254_a.sobj'
name = 'uncultured.thiotrichales.bacterium_AY934254_a'
P = RNAPolytope.construct_from_file(sobj_name)
pl=plot_suboptimal_density_slice(name,P,colormap=my_Greys)

pl.save('plts/tRNA_subopt_slice.pdf')

Error in lines 1-1 Traceback (most recent call last): File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/smc_sagews/sage_server.py", line 995, in execute exec compile(block+'\n', '', 'single') in namespace, locals File "", line 1, in <module> File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/sage/misc/decorators.py", line 473, in wrapper return func(*args, **kwds) File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/sage/plot/graphics.py", line 3212, in save figure.tight_layout() File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/matplotlib/figure.py", line 1756, in tight_layout self.subplots_adjust(**kwargs) File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/matplotlib/figure.py", line 1612, in subplots_adjust self.subplotpars.update(*args, **kwargs) File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/matplotlib/figure.py", line 226, in update raise ValueError('left cannot be >= right') ValueError: left cannot be >= right

# for i,sobj_name in enumerate(sobj_names):
#     if i==5:
#         name = names[i]
#         P = RNAPolytope.construct_from_file(sobj_name)
#         print "sobj_name='{0}'".format(sobj_name)
#         print "name = '{0}'".format(name)
#         pl=plot_suboptimal_density_slice(name,P,colormap=my_Greys)
#         break
#     if i==10:
#         break

# for i,sobj_name in enumerate(sobj_names):
#     if i==4:
#         name = names[i]
#         P = RNAPolytope.construct_from_file(sobj_name)
#         plot_max_avg_acc_slice(name,P,my_Greys).show(figsize=12)
#     if i==10:
#         break

#for i,sobj_name in enumerate(sobj_names):
#    name = names[i]
#    P = RNAPolytope.construct_from_file(sobj_name)
#    plot_suboptimal_density_slice(name,P)
#    if i==5:
#        break

0.01875**3

6.59179687500000e-6

def plot_polyhedron(p,facecolor='green',edgecolor='blue',vertexcolor='red',faces=True,edges=True,vertices=True,faceopacity=0.25,edgeopacity=1,vertexopacity=1):
    vxs = p.faces(0)
    edgs = p.faces(1)
    facs = p.faces(2)
    P = Graphics()
    if faces:
        for f in facs:
            P+=polygon3d(f.vertices(),color=facecolor,alpha=faceopacity)
    if edges:
        for f in edgs:
            P+=line3d(f.vertices(),color=edgecolor,alpha=edgeopacity)
    if vertices:
        for f in vxs:
            P+=point3d(f.vertices(),color=vertexcolor,alpha=vertexopacity)
    return P
def plot_slices(slices):
    P = Graphics()
    slices = sorted(slices,key=lambda sl:sl.volume(),reverse=True)
    n= len(slices)
    max_vol = slices[0].volume()
    for i,sl in enumerate(slices):
        rgb = my_hot(float(i)/n)#Color([random() for i in range(3)])
        P+=plot_polyhedron(sl,edgeopacity=0.5,facecolor=rgb,faceopacity=0.1,vertexopacity=0.3)
    return P

#plot_slices(sl_part[0])+plot_polyhedron(inter,faceopacity=1)

# max_name = None
# max_P = None
# max_subopt=0
# for i,sobj_name in enumerate(sobj_names):
#     name = names[i]
#     P = RNAPolytope.construct_from_file(sobj_name)
#     region_density, totals, maxes = get_suboptimal_ranking(name,P,0)
#     if maxes[0]>max_subopt:
#         max_name = name
#         max_P = P
#         max_subopt = maxes[0]
# plot_suboptimal_density_slice(max_name,max_P)

G

Error in lines 1-1 Traceback (most recent call last): File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/smc_sagews/sage_server.py", line 995, in execute exec compile(block+'\n', '', 'single') in namespace, locals File "", line 1, in <module> NameError: name 'G' is not defined

tweak = sorted(sl_part[0],key=lambda x:x.volume(),reverse=true)

for sl in tweak:#sl_part[0]:
    print "{0}".format(slice2name[sl])
    print "Volume: {0: 8.2f}; center: ({1: 8.2f},{2: 8.2f},{3: 8.2f})".format(float(sl.volume()),float(sl.center()[0]),float(sl.center()[1]),float(sl.center()[2]))
    inter = inter.intersection(sl)

sl = inter

%time
mc=G.clique_maximum()

CPU time: 0.00 s, Wall time: 0.00 s

H=G.copy()
H.delete_vertices(mc)

mc2=H.clique_maximum()

mc2

[12, 17, 24, 30]

I=H.copy()
I.delete_vertices(mc2)
mc3 = I.clique_maximum()

mc3

[7, 15, 16, 28]

inter = l2s[mc[0]]
for i in mc[1:]:
    inter = inter.intersection(l2s[i])

inter

A 3-dimensional polyhedron in QQ^3 defined as the convex hull of 4 vertices

inter.volume()

9/8000

len(mc)

15

inter.vertices_list()

[[129/10, 3/20, -3], [12, 3/16, -57/20], [123/10, 3/20, -3], [117/10, 3/20, -27/10]]

def get_bounding_box(p):
    xmin = min([v[0] for v in p.vertices_list() ])
    xmax = max([v[0] for v in p.vertices_list() ])
    ymin = min([v[1] for v in p.vertices_list() ])
    ymax = max([v[1] for v in p.vertices_list() ])
    zmin = min([v[2] for v in p.vertices_list() ])
    zmax = max([v[2] for v in p.vertices_list() ])
    return ((xmin,ymin,zmin),(xmax,ymax,zmax))
get_bounding_box(inter)

((117/10, 3/20, -3), (129/10, 3/16, -27/10))

p=inter.plot()

print len(G.clique_maximum())
print len(G.complement().clique_maximum())

15 14

14

cl

[[0, 1, 3, 10, 11, 14, 20, 21, 23, 26, 27, 29, 33, 34, 37], [12, 17, 24, 30], [7, 15, 16, 28], [2, 22, 31], [18, 38], [8, 32], [4, 19], [36], [35], [25], [13], [9], [6], [5]]

def non_empty_intersection(slices):
    inter = slices[0]
    for sl in slices[1:]:
        inter = inter.intersection(sl)
        if len(inter.rays())>0 or inter.volume()>0:
            continue
        return False
    return True

slice_part = [[l2s[x] for x in c] for c in cl]

[non_empty_intersection(p) for p in slice_part]

[True, True, True, True, True, True, True, True, True, True, True, True, True, True]

range(10)[0:]

[0, 1, 2, 3, 4, 5, 6, 7, 8, 9]

[n for n in names if n.find('uncultured.thiotri')>-1]

['uncultured.thiotrichales.bacterium_AY934254_a']

def find_diff_suboptimal_structs(d):
#    n='uncultured.thiotrichales.bacterium_AY934254_a'
    for n in names:
        sn = get_sobj_name(n)
        print sn
        P=RNAPolytope.construct_from_file(sn)

        for sl in P.d1_slices():
            subopts,sign = get_subopt_structures(n,sl.original_structure,sign=True)
            sss = set([s.replace('<','.').replace('>','.') for s in subopts])
            if len(sss)>=d and sign==(0,0,0):
                print list(sss)
                print sign
                break
find_diff_suboptimal_structs(4)

sobjs/nitrosomonas.europaea.atcc.19718_BX321863_b.sobj sobjs/lactobacillus.casei.atcc.334_CP000423_j.sobj sobjs/methanosarcina.barkeri.str..fusaro_CP000099_bi.sobj sobjs/uncultured.alpha.proteobacterium_AY933788.sobj sobjs/tropheryma.whipplei.str..twist_AE014184_h.sobj sobjs/uncultured.thiotrichales.bacterium_AY934254_a.sobj sobjs/synechococcus.sp..wh.8102_BX569691.sobj sobjs/yersinia.pestis.kim_AE009952_ae.sobj sobjs/uncultured.alpha.proteobacterium_DQ093144_a.sobj sobjs/carboxydothermus.hydrogenoformans.z-2901_CP000141_af.sobj sobjs/haloferax.volcanii_HALTRK1.sobj sobjs/homo.sapiens_X04779.sobj sobjs/arabidopsis.thaliana_AC006340.sobj sobjs/chlamydophila.pneumoniae.tw-183_AE017159_j.sobj sobjs/thermoplasma.volcanium.gss1_BA000011_b.sobj sobjs/uncultured.gamma.proteobacterium_AY934184_a.sobj sobjs/nocardia.farcinica.ifm.10152_AP006618_az.sobj sobjs/chlorobium.phaeobacteroides.dsm.266_CP000492_ab.sobj sobjs/scenedesmus.obliquus_DQ396875_n.sobj sobjs/pyrococcus.abyssi.ge5_AJ248288.sobj sobjs/euglena.longa_AJ294725_s.sobj sobjs/marchantia.polymorpha_M.CAU.C.Mar.polym.1.sobj sobjs/lactobacillus.delbrueckii.subsp..bulgaricus.atcc.baa-365_CP000412_cj.sobj sobjs/uncultured.alteromonadales.bacterium_AY934393.sobj sobjs/methanospirillum.hungatei.jf-1_CP000254_ad.sobj sobjs/drosophila.melanogaster_DROTRF2.sobj sobjs/pelobacter.carbinolicus.dsm.2380_CP000142_af.sobj sobjs/bordetella.parapertussis_BX640433_e.sobj sobjs/magnetococcus.sp..mc-1_CP000471_ao.sobj

Error in lines 13-13 Traceback (most recent call last): File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/smc_sagews/sage_server.py", line 995, in execute exec compile(block+'\n', '', 'single') in namespace, locals File "", line 1, in <module> File "", line 7, in find_diff_suboptimal_structs File "", line 2, in get_subopt_structures File "", line 3, in load_obj File "/projects/sage/sage-7.6/local/lib/python/pickle.py", line 1384, in load return Unpickler(file).load() File "/projects/sage/sage-7.6/local/lib/python/pickle.py", line 864, in load dispatch[key](self) File "/projects/sage/sage-7.6/local/lib/python/pickle.py", line 964, in load_binfloat self.append(unpack('>d', self.read(8))[0]) File "src/cysignals/signals.pyx", line 252, in cysignals.signals.python_check_interrupt (build/src/cysignals/signals.c:2854) File "src/cysignals/signals.pyx", line 97, in cysignals.signals.sig_raise_exception (build/src/cysignals/signals.c:1303) KeyboardInterrupt

len(set(['a','b','a']))

2

%md
```>
CCAUAGCUCAGCUGGGAGAGCACCUGCUUUGCAAGCAGGGGGUCGGCGGUUCGACCCCGCCUGGCUCCACCAG
(.......(..(...)...((....))...(....).).(.(....).....(....)....(....)..).)
(.......(..(...)...((....))...(....).).(.(....).....(....)....(....)..).)
(.......(..(...)...((....))...(....).).(.(....).....(....)....(....)..).)
(.......(..(...)...((....))...(....).).(.(....).....(....)....(....)..).)```

%md
# 4 Suboptimal structures for a region
Filename: sobjs/uncultured.thiotrichales.bacterium_AY934254_a.sobj
```
['(.......(..(...)...((....))...(....).).(.(....).....(....)....(....)..).)', '(.......(..(...)...((....))...(....).).(.(....).....(....)....(....)..).)', '(.......(..(...)...((....))...(....).).(.(....).....(....)....(....)..).)', '(.......(..(...)...((....))...(....).).(.(....).....(....)....(....)..).)']
(3, 30, 11, (5, 2))
```

4 Suboptimal structures for a region

Filename: sobjs/uncultured.thiotrichales.bacterium_AY934254_a.sobj

['(.......(..(...)...((....))...(....).).(.(....).....(....)....(....)..).)', '(.......(..(...)...((....))...(....).).(.(....).....(....)....(....)..).)', '(.......(..(...)...((....))...(....).).(.(....).....(....)....(....)..).)', '(.......(..(...)...((....))...(....).).(.(....).....(....)....(....)..).)']
(3, 30, 11, (5, 2))

for sl in P.d1_slices():
    subopts,sign = get_subopt_structures(n,sl.original_structure,sign=True)
    sss= set([s.replace('<','.').replace('>','.') for s in subopts])
    if len(subopts)==2 and sign==(0,0,0):
        print [s.replace('<','.').replace('>','.') for s in subopts]
        print sign
        break

['..........(((.(.(((((....))))).).)))..(((((((((((.......))))).)))))).....', '..........(((.(.(((((....))))).).)))..(((((((((((.......))))).)))))).....'] (0, 0, 0)

'..........(((.(.(((((....))))).).)))..(((((((((((.......))))).)))))).....'== '..........(((.(.(((((....))))).).)))..(((((((((((.......))))).)))))).....'

True

print sl.original_structure

(>...<(...)(...)>.<(...)><(...)>....<(>(...)>.<(...)>...<)>....(...)>..<)

1+1

2

# This is a comment
# follows a print statement
print("Printed message")

Printed message

︠c2b7d1c0-9842-4189-98d2-3c0dda3b9d34︠
D = load_obj("st_tuples/"+name+"_ss.pkl")
K=D.keys()
for k in K:
    if len(D[k])==4:
        S,(x,y,z,w,t1,t2,t3,t4,acc) in D[k]
        print [s[0] for s in D[k]],(x,y,z,w)

sobjs/uncultured.thiotrichales.bacterium_AY934254_a.sobj
['<(>((((...))))((>(>((....))>(((....)))((((((((....))))))))><)>(....)>))>)', '<(><(((...)))(((>((((....))))((....))>((((((((....))))))))><)>(....)>))>)', '<(><(((...)))(((>((((....))))((....))>((((((((....))))))))>))>(....)><)>)']
(3, 9, 10, (-49, 2))