Update bounds with calls from CNVkit, inferred copy numbers and p-values from THetA.
769635
def _update_with_calls(result_file, cnv_file):
    
    results = {}
    with open(result_file) as in_handle:
        in_handle.readline()  # header
        _, _, cs, ps = in_handle.readline().strip().split()
        for i, (c, p) in enumerate(zip(cs.split(":"), ps.split(","))):
            results[i] = (c, p)
    cnvs = {}
    with open(cnv_file) as in_handle:
        for line in in_handle:
            chrom, start, end, _, count = line.rstrip().split()[:5]
            cnvs[(chrom, start, end)] = count
    def update(i, line):
        parts = line.rstrip().split("\t")
        chrom, start, end = parts[1:4]
        parts += cnvs.get((chrom, start, end), ".")
        parts += list(results[i])
        return "\t".join(parts) + "\n"
    return update
734628
def update_terminal_regions(self, tree, X, y, residual, y_pred,
                                sample_weight, sample_mask,
                                learning_rate=1.0, k=0):
        
        # update predictions
        y_pred[:, k] += learning_rate * tree.predict(X).ravel()
351058
def _cnvkit_metrics(cnns, target_bed, antitarget_bed, cov_interval, items):
    
    if cov_interval == "genome":
        return cnns
    target_cnn = [x["file"] for x in cnns if x["cnntype"] == "target"][0]
    background_file = "%s-flatbackground.cnn" % utils.splitext_plus(target_cnn)[0]
    background_file = cnvkit_background([], background_file, items, target_bed, antitarget_bed)
    cnr_file, data = _cnvkit_fix_base(cnns, background_file, items, "-flatbackground")
    cns_file = _cnvkit_segment(cnr_file, cov_interval, data)
    metrics_file = "%s-metrics.txt" % utils.splitext_plus(target_cnn)[0]
    if not utils.file_exists(metrics_file):
        with file_transaction(data, metrics_file) as tx_metrics_file:
            cmd = [_get_cmd(), "metrics", "-o", tx_metrics_file, "-s", cns_file, "--", cnr_file]
            do.run(_prep_cmd(cmd, tx_metrics_file), "CNVkit metrics")
    metrics = _read_metrics_file(metrics_file)
    out = []
    for cnn in cnns:
        cnn["metrics"] = metrics
        out.append(cnn)
    return out
965861
def _cnv_prioritize(data):
    
    supported = {"cnvkit": {"inputs": ["call_file", "segmetrics"], "fn": _cnvkit_prioritize}}
    pcall = None
    priority_files = None
    for call in data.get("sv", []):
        if call["variantcaller"] in supported:
            priority_files = [call.get(x) for x in supported[call["variantcaller"]]["inputs"]]
            priority_files = [x for x in priority_files if x is not None and utils.file_exists(x)]
            if len(priority_files) == len(supported[call["variantcaller"]]["inputs"]):
                pcall = call
                break
    prioritize_by = tz.get_in(["config", "algorithm", "svprioritize"], data)
    if pcall and prioritize_by:
        out_file = "%s-prioritize.tsv" % utils.splitext_plus(priority_files[0])[0]
        gene_list = _find_gene_list_from_bed(prioritize_by, out_file, data)
        if gene_list:
            with open(gene_list) as in_handle:
                genes = [x.strip() for x in in_handle]
            args = [dd.get_sample_name(data), genes] + priority_files
            df = supported[pcall["variantcaller"]]["fn"](*args)
            with file_transaction(data, out_file) as tx_out_file:
                df.to_csv(tx_out_file, sep="\t", index=False)
            pcall["priority"] = out_file
    return data
920996
def _run_cnvkit_shared_orig(inputs, backgrounds):
    
    work_dir = _sv_workdir(inputs[0])
    raw_work_dir = utils.safe_makedir(os.path.join(work_dir, "raw"))
    background_name = dd.get_sample_name(backgrounds[0]) if backgrounds else "flat"
    background_cnn = os.path.join(raw_work_dir, "%s_background.cnn" % (background_name))
    ckouts = []
    for cur_input in inputs:
        cur_raw_work_dir = utils.safe_makedir(os.path.join(_sv_workdir(cur_input), "raw"))
        out_base, out_base_old = _bam_to_outbase(dd.get_align_bam(cur_input), cur_raw_work_dir, cur_input)
        if utils.file_exists(out_base_old + ".cns"):
            out_base = out_base_old
        ckouts.append({"cnr": "%s.cnr" % out_base,
                       "cns": "%s.cns" % out_base})
    if not utils.file_exists(ckouts[0]["cns"]):
        cov_interval = dd.get_coverage_interval(inputs[0])
        samples_to_run = list(zip(["background"] * len(backgrounds), backgrounds)) + \
                         list(zip(["evaluate"] * len(inputs), inputs))
        # New style shared SV bins
        if tz.get_in(["depth", "bins", "target"], inputs[0]):
            target_bed = tz.get_in(["depth", "bins", "target"], inputs[0])
            antitarget_bed = tz.get_in(["depth", "bins", "antitarget"], inputs[0])
            raw_coverage_cnns = reduce(operator.add,
                                       [_get_general_coverage(cdata, itype) for itype, cdata in samples_to_run])
        # Back compatible with pre-existing runs
        else:
            target_bed, antitarget_bed = _get_original_targets(inputs[0])
            raw_coverage_cnns = reduce(operator.add,
                                       [_get_original_coverage(cdata, itype) for itype, cdata in samples_to_run])
        # Currently metrics not calculated due to speed and needing re-evaluation
        # We could re-enable with larger truth sets to evaluate background noise
        # But want to reimplement in a more general fashion as part of normalization
        if False:
            coverage_cnns = reduce(operator.add,
                                [_cnvkit_metrics(cnns, target_bed, antitarget_bed, cov_interval,
                                                    inputs + backgrounds)
                                    for cnns in tz.groupby("bam", raw_coverage_cnns).values()])
            background_cnn = cnvkit_background(_select_background_cnns(coverage_cnns),
                                                background_cnn, inputs, target_bed, antitarget_bed)
        else:
            coverage_cnns = raw_coverage_cnns
            background_cnn = cnvkit_background([x["file"] for x in coverage_cnns if x["itype"] == "background"],
                                                background_cnn, inputs, target_bed, antitarget_bed)
        parallel = {"type": "local", "cores": dd.get_cores(inputs[0]), "progs": ["cnvkit"]}
        fixed_cnrs = run_multicore(_cnvkit_fix,
                                   [(cnns, background_cnn, inputs, ckouts) for cnns in
                                    tz.groupby("bam", [x for x in coverage_cnns
                                                       if x["itype"] == "evaluate"]).values()],
                                   inputs[0]["config"], parallel)
        [_cnvkit_segment(cnr, cov_interval, data, inputs + backgrounds) for cnr, data in fixed_cnrs]
    return ckouts
789098
def update_bounds(self, bounds):
        '''

        '''
        starts =  bounds[:,0,:]
        ends = bounds[:,1,:]

        self.bounds = bounds
        self.lengths =  np.sqrt(((ends - starts)**2).sum(axis=1))
        
        vertices, normals, colors = self._process_reference()
        
        self.tr.update_vertices(vertices)
        self.tr.update_normals(normals)
1122254
def update_plot_limits(ax, white_space):
    

    if hasattr(ax, 'zz_dataLim'):
        bounds = ax.xy_dataLim.bounds
        ax.set_xlim(bounds[0] - white_space, bounds[0] + bounds[2] + white_space)
        ax.set_ylim(bounds[1] - white_space, bounds[1] + bounds[3] + white_space)

        bounds = ax.zz_dataLim.bounds
        ax.set_zlim(bounds[0] - white_space, bounds[0] + bounds[2] + white_space)
    else:
        bounds = ax.dataLim.bounds
        assert not any(map(np.isinf, bounds)), 'Cannot set bounds if dataLim has infinite elements'
        ax.set_xlim(bounds[0] - white_space, bounds[0] + bounds[2] + white_space)
        ax.set_ylim(bounds[1] - white_space, bounds[1] + bounds[3] + white_space)
1065473
def _get_calls(data, cnv_only=False):
    
    cnvs_supported = set(["cnvkit", "battenberg"])
    out = {}
    for sv in data.get("sv", []):
        if not cnv_only or sv["variantcaller"] in cnvs_supported:
            out[sv["variantcaller"]] = sv
    return out
1167800
def _param_bounds(self, gamma, q):
        
        if gamma < 1.4:
            gamma = 1.4
        if gamma > 2.9:
            gamma = 2.9
        if q < 0.01:
            q = 0.01
        return float(gamma), q
1261365
private void updateThreshold(int i)
    {
        double Fi = fcache[i];
        double F_tilde_i = b_low;
        
        if (I0_b[i] || I2[i])
            F_tilde_i = Fi + epsilon;
        else if (I0_a[i] || I1[i])
            F_tilde_i = Fi - epsilon;
        
        double F_bar_i = b_up;
        if (I0_a[i] || I3[i])
            F_bar_i = Fi - epsilon;
        else if (I0_b[i] || I1[i])
            F_bar_i = Fi + epsilon;
        
        //update the bounds
        
        if (b_low < F_tilde_i)
        {
            b_low = F_tilde_i;
            i_low = i;
        }
        
        if (b_up > F_bar_i)
        {
            b_up = F_bar_i;
            i_up = i;
        }
    }
688485
def p_funcvardecls(self, p):
        ''
        p[0] = p[1] + (p[2],)
        p.set_lineno(0, p.lineno(1))
127115
function fp_inv(r, n, p, l) {
	if (mpn_zero_p(n, l)) {
		throw new RangeError("not invertible");
	}
	var u = _fp_u, v = _fp_v, s = _fp_s;
	mpn_copyi(u, n, l), mpn_copyi(v, p, l);
	mpn_zero(r, l), mpn_zero(s, l);
	r[0] = 1;
	for (;;) {
		if (mpn_one_p(u, l)) {
			return r;
		}
		if (mpn_one_p(v, l)) {
			mpn_copyi(r, s, l);
			return r;
		}
		while (mpn_even_p(u, l)) {
			mpn_shr1(u, u, l);
			if (mpn_even_p(r, l)) {
				mpn_shr1(r, r, l);
			}
			else {
				var c = mpn_add(r, r, p, l) << 30;
				mpn_shr1(r, r, l);
				r[l - 1] |= c;
			}
		}
		while (mpn_even_p(v, l)) {
			mpn_shr1(v, v, l);
			if (mpn_even_p(s, l)) {
				mpn_shr1(s, s, l);
			}
			else {
				var c = mpn_add(s, s, p, l) << 30;
				mpn_shr1(s, s, l);
				s[l - 1] |= c;
			}
		}
		if (mpn_cmp(u, v, l) >= 0) {
			mpn_sub(u, u, v, l);
			fp_sub(r, r, s, p, l);
		}
		else {
			mpn_sub(v, v, u, l);
			fp_sub(s, s, r, p, l);
		}
	}
}
1084398
def optimize(self, x0, f=None, df=None, f_df=None):
        

        if len(self.bounds) == 1:
            raise IndexError("CMA does not work in problems of dimension 1.")

        try:
            import cma
            def CMA_f_wrapper(f):
                def g(x):
                    return f(np.array([x]))
                return g
            lB = np.asarray(self.bounds)[:,0]
            uB = np.asarray(self.bounds)[:,1]
            x = cma.fmin(CMA_f_wrapper(f), x0, 0.6, options={"bounds":[lB, uB], "verbose":-1})[0]
            return np.atleast_2d(x), f(np.atleast_2d(x))
        except ImportError:
            print("Cannot find cma library, please install it to use this option.")
            raise
210774
def _run_theta(cnv_info, data, work_dir, run_n3=True):
    
    out = {"caller": "theta"}
    max_normal = "0.9"
    opts = ["-m", max_normal]
    n2_result = _safe_run_theta(cnv_info["theta_input"], os.path.join(work_dir, "n2"), ".n2.results",
                                ["-n", "2"] + opts, data)
    if n2_result:
        out["estimate"] = n2_result
        if run_n3:
            n2_bounds = "%s.withBounds" % os.path.splitext(n2_result)[0]
            n3_result = _safe_run_theta(n2_bounds, os.path.join(work_dir, "n3"), ".n3.results",
                                        ["-n", "3", "--RESULTS", n2_result] + opts,
                                        data)
            if n3_result:
                best_result = _select_model(n2_bounds, n2_result, n3_result,
                                            os.path.join(work_dir, "n3"), data)
                out["estimate"] = best_result
                out["cnvs"] = _merge_theta_calls(n2_bounds, best_result, cnv_info["vrn_file"], data)
    return out
1254542
def _update_3d_datalim(ax, obj):
    ''''''
    min_bounding_box, max_bounding_box = geom.bounding_box(obj)
    xy_bounds = np.vstack((min_bounding_box[:COLS.Z],
                           max_bounding_box[:COLS.Z]))
    ax.xy_dataLim.update_from_data_xy(xy_bounds, ignore=False)

    z_bounds = np.vstack(((min_bounding_box[COLS.Z], min_bounding_box[COLS.Z]),
                          (max_bounding_box[COLS.Z], max_bounding_box[COLS.Z])))
    ax.zz_dataLim.update_from_data_xy(z_bounds, ignore=False)
342186
def _associate_cnvkit_out(ckouts, items, is_somatic=False):
    
    assert len(ckouts) == len(items)
    out = []
    upload_counts = collections.defaultdict(int)
    for ckout, data in zip(ckouts, items):
        ckout = copy.deepcopy(ckout)
        ckout["variantcaller"] = "cnvkit"
        if utils.file_exists(ckout["cns"]) and _cna_has_values(ckout["cns"]):
            ckout = _add_seg_to_output(ckout, data)
            ckout = _add_gainloss_to_output(ckout, data)
            ckout = _add_segmetrics_to_output(ckout, data)
            ckout = _add_variantcalls_to_output(ckout, data, items, is_somatic)
            # ckout = _add_coverage_bedgraph_to_output(ckout, data)
            ckout = _add_cnr_bedgraph_and_bed_to_output(ckout, data)
            if "svplots" in dd.get_tools_on(data):
                ckout = _add_plots_to_output(ckout, data)
            ckout["do_upload"] = upload_counts[ckout.get("vrn_file")] == 0
        if "sv" not in data:
            data["sv"] = []
        data["sv"].append(ckout)
        if ckout.get("vrn_file"):
            upload_counts[ckout["vrn_file"]] += 1
        out.append(data)
    return out
833702
def _bnd(self, xloc, dist, cache):
        """"""
        return numpy.log(evaluation.evaluate_bound(
            dist, numpy.e**xloc, cache=cache))