Python Euclid.adjacency_to_laplacian方法代码示例

# 需要导入模块: import Euclid [as 别名]
# 或者: from Euclid import adjacency_to_laplacian [as 别名]
def get_response_content(fs):
    locations = get_locations()
    np_locs = [np.array(p) for p in locations]
    edges = get_edges()
    npoints = len(locations)
    # start writing the response
    np.set_printoptions(linewidth=200)
    out = StringIO()
    # print the layout data
    print >> out, 'POINTS'
    for i, (x, y) in enumerate(locations):
        print >> out, i, x, y
    print >> out, 'EDGES'
    for i, j in edges:
        print >> out, i, j
    print >> out
    # show the unweighted adjacency matrix
    UA = np.zeros((npoints, npoints))
    for i, j in edges:
        UA[i, j] = 1
        UA[j, i] = 1
    print >> out, 'unweighted adjacency matrix:'
    print >> out, UA
    print >> out
    # show the unweighted laplacian matrix
    UL = Euclid.adjacency_to_laplacian(UA)
    print >> out, 'unweighted laplacian matrix:'
    print >> out, UL
    print >> out
    # show the weighted adjacency matrix
    WA = np.zeros((npoints, npoints))
    for i, j in edges:
        d = np.linalg.norm(np_locs[i] - np_locs[j]) / math.sqrt(2.0)
        w = 1.0 / d
        WA[i, j] = w
        WA[j, i] = w
    print >> out, 'weighted adjacency matrix:'
    print >> out, WA
    print >> out
    # show the weighted laplacian matrix
    WL = Euclid.adjacency_to_laplacian(WA)
    print >> out, 'weighted laplacian matrix:'
    print >> out, WL
    print >> out
    # remove the two internal nodes by schur complementation
    ntips = 4
    schur_L = SchurAlgebra.schur_helper(WL, 2)
    X = Euclid.dccov_to_points(np.linalg.pinv(schur_L))
    print >> out, 'schur graph layout:'
    print >> out, 'POINTS'
    for i, v in enumerate(X):
        print >> out, i, v[0], v[1]
    print >> out, 'EDGES'
    for i in range(ntips):
        for j in range(i+1, ntips):
            print >> out, i, j
    # return the response
    return out.getvalue()

# 需要导入模块: import Euclid [as 别名]
# 或者: from Euclid import adjacency_to_laplacian [as 别名]
def get_response_content(fs):
    # build the newick tree from the string
    tree = NewickIO.parse(fs.tree_string, FelTree.NewickTree)
    nvertices = len(list(tree.preorder()))
    nleaves = len(list(tree.gen_tips()))
    # get ordered ids with the leaves first
    ordered_ids = get_ordered_ids(tree)
    # get the adjacency matrix and the augmented adjacency matrix
    A = np.array(tree.get_affinity_matrix(ordered_ids))
    A_aug = get_augmented_adjacency(A, nleaves, fs.ndups, fs.strength)
    # get the laplacian matrices
    L = Euclid.adjacency_to_laplacian(A)
    L_aug = Euclid.adjacency_to_laplacian(A_aug)
    # get the schur complement
    R = SchurAlgebra.mschur(L, set(range(nleaves, nvertices)))
    R_pinv = np.linalg.pinv(R)
    vals, vecs = EigUtil.eigh(R_pinv)
    # get the scaled Fiedler vector for the Schur complement
    w, v = EigUtil.principal_eigh(R_pinv)
    fiedler = v * math.sqrt(w)
    # get the eigendecomposition of the augmented Laplacian
    L_aug_pinv = np.linalg.pinv(L_aug)
    vals_aug, vecs_aug = EigUtil.eigh(L_aug_pinv)
    # get the scaled Fiedler vector for the augmented Laplacian
    w_aug, v_aug = EigUtil.principal_eigh(L_aug_pinv)
    fiedler_aug = v_aug * math.sqrt(w_aug)
    # report the results
    np.set_printoptions(linewidth=300)
    out = StringIO()
    print >> out, 'Laplacian matrix:'
    print >> out, L
    print >> out
    print >> out, 'Schur complement of Laplacian matrix:'
    print >> out, R
    print >> out
    print >> out, 'scaled Fiedler vector of Schur complement:'
    print >> out, fiedler
    print >> out
    print >> out, 'eigenvalues of pinv of Schur complement:'
    print >> out, vals
    print >> out
    print >> out, 'corresponding eigenvectors of pinv of Schur complement:'
    print >> out, np.array(vecs).T
    print >> out
    print >> out
    print >> out, 'augmented Laplacian matrix:'
    print >> out, L_aug
    print >> out
    print >> out, 'scaled Fiedler vector of augmented Laplacian:'
    print >> out, fiedler_aug
    print >> out
    print >> out, 'eigenvalues of pinv of augmented Laplacian:'
    print >> out, vals_aug
    print >> out
    print >> out, 'rows are eigenvectors of pinv of augmented Laplacian:'
    print >> out, np.array(vecs_aug)
    return out.getvalue()

# 需要导入模块: import Euclid [as 别名]
# 或者: from Euclid import adjacency_to_laplacian [as 别名]
def create_laplacian_matrix(lena, lenb, lenc):
    """
    @param lena: integer length of first branch.
    @param lenb: integer length of second branch.
    @param lenc: integer length of third branch.
    """
    N = 1 + lena + lenb + lenc
    A = np.zeros((N,N), dtype=float)
    # Add connections to the hub vertex.
    if lena:
        A[0, 1] = 1
        A[1, 0] = 1
    if lenb:
        A[0, lena+1] = 1
        A[lena+1, 0] = 1
    if lenc:
        A[0, lena+lenb+1] = 1
        A[lena+lenb+1, 0] = 1
    # Add tridiagonal connections on the first branch.
    for i in range(lena-1):
        j = i + 1
        A[j, j+1] = 1
        A[j+1, j] = 1
    # Add tridiagonal connections on the second branch.
    for i in range(lenb-1):
        j = lena + i + 1
        A[j, j+1] = 1
        A[j+1, j] = 1
    # Add tridiagonal connections on the second branch.
    for i in range(lenc-1):
        j = lena + lenb + i + 1
        A[j, j+1] = 1
        A[j+1, j] = 1
    L = Euclid.adjacency_to_laplacian(A)
    return L

# 需要导入模块: import Euclid [as 别名]
# 或者: from Euclid import adjacency_to_laplacian [as 别名]
def get_full_tree_message(tree, m_to_string):
    """
    In this function we find the Fiedler split of the full tree.
    @param tree: each node in this tree must have a name
    @param m_to_string: a function that converts a matrix to a string
    @return: a message about the split of the tips of the tree induced by the fiedler vector
    """
    out = StringIO()
    # get the alphabetically ordered names
    ordered_names = list(sorted(node.get_name() for node in tree.preorder()))
    # get the corresponding ordered ids
    name_to_id = dict((node.get_name(), id(node)) for node in tree.preorder())
    ordered_ids = [name_to_id[name] for name in ordered_names]
    # get the full weighted adjacency matrix
    A = np.array(tree.get_affinity_matrix(ordered_ids))
    print >> out, 'the weighted reciprocal adjacency matrix of the full tree:'
    print >> out, m_to_string(get_reciprocal_matrix(A))
    print >> out
    # get the full Laplacian matrix
    L = Euclid.adjacency_to_laplacian(A)
    # get the fiedler split
    v = BuildTreeTopology.laplacian_to_fiedler(L)
    print >> out, 'the Fiedler split of the full tree:'
    for name, value in zip(ordered_names, v):
        print >> out, name, ':', value
    return out.getvalue().strip()

# 需要导入模块: import Euclid [as 别名]
# 或者: from Euclid import adjacency_to_laplacian [as 别名]
def get_response_content(fs):
    out = StringIO()
    # try to make some graphs
    unconnected_count = 0
    invalid_split_count = 0
    valid_split_count = 0
    for graph_index in range(fs.ngraphs):
        G = erdos_renyi(fs.nvertices, fs.pedge)
        if is_connected(G):
            # add interesting edge weights
            add_exponential_weights(G)
            # turn the adjacency matrix into a laplacian matrix
            L = Euclid.adjacency_to_laplacian(G)
            for v in range(fs.nvertices):
                small_index_to_big_index = {}
                for i_small, i_big in enumerate([i for i in range(fs.nvertices) if i != v]):
                    small_index_to_big_index[i_small] = i_big
                # take the schur complement with respect to the given vertex
                L_reduced = get_single_element_schur_complement(L, v)
                assert len(L_reduced) == len(L) - 1
                # get the loadings of the vertices of the reduced graph
                if fs.fiedler_cut:
                    Y_reduced = BuildTreeTopology.laplacian_to_fiedler(L_reduced)
                elif fs.random_cut:
                    Y_reduced = get_random_vector(L_reduced)
                assert len(Y_reduced) == len(L_reduced)
                # expand the fiedler vector with positive and negative valuations for the removed vertex
                found_valid_split = False
                for augmented_loading in (-1.0, 1.0):
                    # get the augmented split vector for this assignment of the removed vertex
                    Y_full = [0]*len(G)
                    for i_reduced, loading in enumerate(Y_reduced):
                        i_big = small_index_to_big_index[i_reduced]
                        Y_full[i_big] = loading
                    Y_full[v] = augmented_loading
                    assert len(Y_full) == len(G)
                    # get the two graphs defined by the split
                    subgraph_a, subgraph_b = list(gen_subgraphs(G, Y_full))
                    # if the subgraphs are both connected then the split is valid
                    if is_connected(subgraph_a) and is_connected(subgraph_b):
                        found_valid_split = True
                # if a valid split was not found then show the matrix
                if found_valid_split:
                    valid_split_count += 1
                else:
                    print >> out, 'Found a matrix that was split incompatibly by a cut of its schur complement!'
                    print >> out, 'matrix:'
                    print >> out, MatrixUtil.m_to_string(G)
                    print >> out, 'index that was removed:', v
                    invalid_split_count += 1
        else:
            unconnected_count += 1
    # show the number of connected and of unconnected graphs
    print >> out, 'this many random graphs were connected:', fs.ngraphs - unconnected_count
    print >> out, 'this many random graphs were not connected:', unconnected_count
    print >> out, 'this many splits were valid:', valid_split_count
    print >> out, 'this many splits were invalid:', invalid_split_count
    # return the result
    return out.getvalue()

# 需要导入模块: import Euclid [as 别名]
# 或者: from Euclid import adjacency_to_laplacian [as 别名]
def get_response_content(fs):
    # get the tree
    tree = NewickIO.parse(fs.tree_string, FelTree.NewickTree)
    # get information about the tree topology
    internal = [id(node) for node in tree.gen_internal_nodes()]
    tips = [id(node) for node in tree.gen_tips()]
    vertices = internal + tips
    ntips = len(tips)
    ninternal = len(internal)
    nvertices = len(vertices)
    # get the ordered ids with the leaves first
    ordered_ids = vertices
    # get the full weighted adjacency matrix
    A = np.array(tree.get_affinity_matrix(ordered_ids))
    # compute the weighted adjacency matrix of the decorated tree
    p = ninternal
    q = ntips
    N = fs.N
    if fs.weight_n:
        weight = float(N)
    elif fs.weight_sqrt_n:
        weight = math.sqrt(N)
    A_aug = get_A_aug(A, weight, p, q, N)
    # compute the weighted Laplacian matrix of the decorated tree
    L_aug = Euclid.adjacency_to_laplacian(A_aug)
    # compute the eigendecomposition
    w, vt = np.linalg.eigh(L_aug)
    # show the output
    np.set_printoptions(linewidth=1000, threshold=10000)
    out = StringIO()
    if fs.lap:
        print >> out, 'Laplacian of the decorated tree:'
        print >> out, L_aug
        print >> out
    if fs.eigvals:
        print >> out, 'eigenvalues:'
        for x in w:
            print >> out, x
        print >> out
    if fs.eigvecs:
        print >> out, 'eigenvector matrix:'
        print >> out, vt
        print >> out
    if fs.compare:
        # get the distance matrix for only the original tips
        D_tips = np.array(tree.get_partial_distance_matrix(tips))
        X_tips = Euclid.edm_to_points(D_tips)
        # wring the approximate points out of the augmented tree
        X_approx = vt[p:p+q].T[1:1+q-1].T / np.sqrt(w[1:1+q-1])
        # do the comparison
        print >> out, 'points from tip-only MDS:'
        print >> out, X_tips
        print >> out
        print >> out, 'approximate points from decorated tree:'
        print >> out, X_approx
        print >> out
    return out.getvalue()

# 需要导入模块: import Euclid [as 别名]
# 或者: from Euclid import adjacency_to_laplacian [as 别名]
def get_response_content(fs):
    # build the newick tree from the string
    tree = NewickIO.parse(fs.tree_string, FelTree.NewickTree)
    nvertices = len(list(tree.preorder()))
    nleaves = len(list(tree.gen_tips()))
    # get ordered ids with the leaves first
    ordered_ids = get_ordered_ids(tree)
    # get the adjacency matrix and the augmented adjacency matrix
    A = np.array(tree.get_affinity_matrix(ordered_ids))
    A_aug = get_augmented_adjacency(A, nleaves, fs.strength)
    # get the laplacian matrices
    L = Euclid.adjacency_to_laplacian(A)
    L_aug = Euclid.adjacency_to_laplacian(A_aug)
    # get the schur complements
    R = SchurAlgebra.mschur(L, set(range(nleaves, nvertices)))
    R_aug = SchurAlgebra.mschur(L_aug, set(range(nleaves, nvertices)))
    # get the scaled Fiedler vectors
    w, v = EigUtil.principal_eigh(np.linalg.pinv(R))
    fiedler = v * math.sqrt(w)
    w_aug, v_aug = EigUtil.principal_eigh(np.linalg.pinv(R_aug))
    fiedler_aug = v_aug * math.sqrt(w_aug)
    # report the results
    np.set_printoptions(linewidth=200)
    out = StringIO()
    print >> out, "Laplacian matrix:"
    print >> out, L
    print >> out
    print >> out, "Schur complement of Laplacian matrix:"
    print >> out, R
    print >> out
    print >> out, "scaled Fiedler vector:"
    print >> out, fiedler
    print >> out
    print >> out, "augmented Laplacian matrix:"
    print >> out, L_aug
    print >> out
    print >> out, "Schur complement of augmented Laplacian matrix:"
    print >> out, R_aug
    print >> out
    print >> out, "scaled Fiedler vector of augmented matrix:"
    print >> out, fiedler_aug
    print >> out
    return out.getvalue()

# 需要导入模块: import Euclid [as 别名]
# 或者: from Euclid import adjacency_to_laplacian [as 别名]
def get_response_content(fs):
    # use a fixed seed if requested
    if fs.seed:
        random.seed(fs.seed)
    # define the max number of rejection iterations
    limit = fs.npoints * 100
    # validate input
    if fs.axis < 0:
        raise ValueError("the mds axis must be nonnegative")
    # get points defining the boundary of africa
    nafrica = len(g_africa_poly)
    africa_edges = [(i, (i + 1) % nafrica) for i in range(nafrica)]
    # get some points and edges inside africa
    points = sample_with_rejection(fs.npoints, g_africa_poly, limit)
    x_list, y_list = zip(*points)
    tri = Triangulation(x_list, y_list)
    tri_edges = [(i + nafrica, j + nafrica) for i, j in tri.edge_db.tolist()]
    # get the whole list of points
    allpoints = g_africa_poly + points
    # refine the list of edges
    tri_edges = list(gen_noncrossing_edges(tri_edges, africa_edges, allpoints))
    tri_edges = get_mst(tri_edges, allpoints)
    alledges = africa_edges + tri_edges
    # make the graph laplacian
    A = np.zeros((len(points), len(points)))
    for ia, ib in tri_edges:
        xa, ya = allpoints[ia]
        xb, yb = allpoints[ib]
        d = math.hypot(xb - xa, yb - ya)
        A[ia - nafrica, ib - nafrica] = 1 / d
        A[ib - nafrica, ia - nafrica] = 1 / d
    L = Euclid.adjacency_to_laplacian(A)
    ws, vs = EigUtil.eigh(np.linalg.pinv(L))
    if fs.axis >= len(ws):
        raise ValueError("choose a smaller mds axis")
    v = vs[fs.axis]
    # get the color and sizes for the points
    v /= max(np.abs(v))
    colors = [(0, 0, 0)] * nafrica + [get_color(x) for x in v]
    radii = [2] * nafrica + [5 for p in points]
    # get the width and height of the drawable area of the image
    width = fs.total_width - 2 * fs.border
    height = fs.total_height - 2 * fs.border
    if width < 1 or height < 1:
        msg = "the image dimensions do not allow for enough drawable area"
        raise HandlingError(msg)
    # draw the image
    ext = Form.g_imageformat_to_ext[fs.imageformat]
    try:
        helper = ImgHelper(allpoints, alledges, fs.total_width, fs.total_height, fs.border)
        return helper.get_image_string(colors, radii, ext)
    except CairoUtil.CairoUtilError as e:
        raise HandlingError(e)

# 需要导入模块: import Euclid [as 别名]
# 或者: from Euclid import adjacency_to_laplacian [as 别名]
def create_laplacian_matrix(nvertices):
    """
    @param affinity: affinity between adjacent vertices
    @param nvertices: the number of vertices in the graph
    @return: a numpy matrix
    """
    affinity = nvertices * 2.0
    A = np.zeros((nvertices, nvertices), dtype=float)
    for i, j in iterutils.pairwise(range(nvertices)):
        A[i,j] = affinity
        A[j,i] = affinity
    L = Euclid.adjacency_to_laplacian(A)
    return L

# 需要导入模块: import Euclid [as 别名]
# 或者: from Euclid import adjacency_to_laplacian [as 别名]
def get_form():
    """
    @return: the body of a form
    """
    # define the default matrix string
    A = np.array(StoerWagner.g_stoer_wagner_affinity)
    M = -Euclid.adjacency_to_laplacian(A)
    # define the form objects
    form_objects = [
            Form.Matrix('matrix', 'matrix', M),
            Form.CheckGroup('drawing_options', 'drawing options', [
                Form.CheckItem('axes', 'draw axes', True),
                Form.CheckItem('connections', 'draw connections', True),
                Form.CheckItem('vertices', 'draw vertices', True)]),
            Form.ImageFormat()]
    return form_objects

# 需要导入模块: import Euclid [as 别名]
# 或者: from Euclid import adjacency_to_laplacian [as 别名]
 def do_search(self, nseconds, sampling_function):
     """
     @param nseconds: allowed search time or None
     @param sampling_function: a function that samples a branch length
     @return: True if a tree was found that met the criteria
     """
     if not self.is_initialized():
         raise RuntimeError("the search was not sufficiently initialized")
     true_splits = self.tree.get_nontrivial_splits()
     start_time = time.time()
     while True:
         elapsed_time = time.time() - start_time
         if nseconds and elapsed_time > nseconds:
             return False
         # assign new sampled branch lengths
         for branch in self.tree.get_branches():
             branch.length = sampling_function()
         # get the distance matrix so we can use a library function to get the split
         D = np.array(self.tree.get_distance_matrix())
         ntips = len(D)
         # get the Laplacian matrix of the full tree and the corresponding Fiedler split of the leaves
         if self.force_difference or self.informative_full_split:
             A_aug = np.array(self.tree.get_weighted_adjacency_matrix(self.id_to_index))
             L_aug = Euclid.adjacency_to_laplacian(A_aug)
             v_aug = BuildTreeTopology.laplacian_to_fiedler(L_aug)
             left_aug, right_aug = BuildTreeTopology.eigenvector_to_split(v_aug)
             left = [x for x in left_aug if x in range(ntips)]
             right = [x for x in right_aug if x in range(ntips)]
             leaf_eigensplit_aug = BuildTreeTopology.make_split(left, right)
             if self.force_difference:
                 if leaf_eigensplit_aug == self.desired_primary_split:
                     self.aug_split_collision_count += 1
                     continue
             if self.informative_full_split:
                 if min(len(s) for s in leaf_eigensplit_aug) < 2:
                     self.aug_split_degenerate_count += 1
                     continue
         # get the eigensplit
         try:
             eigensplit = BuildTreeTopology.split_using_eigenvector(D)
         except BuildTreeTopology.DegenerateSplitException, e:
             self.degenerate_primary_split_count += 1
             continue
         except BuildTreeTopology.InvalidSpectralSplitException, e:
             self.error_primary_split_count += 1
             continue

# 需要导入模块: import Euclid [as 别名]
# 或者: from Euclid import adjacency_to_laplacian [as 别名]
def main(fs):
    # use a fixed seed if requested
    if fs.seed:
        random.seed(fs.seed)
    # define the max number of rejection iterations
    limit = fs.npoints * 100
    # validate input
    if fs.axis < 0:
        raise ValueError("the mds axis must be nonnegative")
    # get points defining the boundary of africa
    nafrica = len(g_africa_poly)
    africa_edges = [(i, (i + 1) % nafrica) for i in range(nafrica)]
    # get some points and edges inside africa
    points = sample_with_rejection(fs.npoints, g_africa_poly, limit)
    x_list, y_list = zip(*points)
    tri = Triangulation(x_list, y_list)
    tri_edges = [(i + nafrica, j + nafrica) for i, j in tri.edge_db.tolist()]
    # get the whole list of points
    allpoints = g_africa_poly + points
    # refine the list of edges
    tri_edges = list(gen_noncrossing_edges(tri_edges, africa_edges, allpoints))
    tri_edges = get_mst(tri_edges, allpoints)
    alledges = africa_edges + tri_edges
    # make the graph laplacian
    A = np.zeros((len(points), len(points)))
    for ia, ib in tri_edges:
        xa, ya = allpoints[ia]
        xb, yb = allpoints[ib]
        d = math.hypot(xb - xa, yb - ya)
        A[ia - nafrica, ib - nafrica] = 1 / d
        A[ib - nafrica, ia - nafrica] = 1 / d
    L = Euclid.adjacency_to_laplacian(A)
    ws, vs = EigUtil.eigh(np.linalg.pinv(L))
    if fs.axis >= len(ws):
        raise ValueError("choose a smaller mds axis")
    v = vs[fs.axis]
    # get the color and sizes for the points
    v /= max(np.abs(v))
    # draw the picture
    helper = ImgHelper(allpoints, alledges, fs.total_width, fs.total_height, fs.border)
    helper.draw_contour_plot(v, nafrica)

# 需要导入模块: import Euclid [as 别名]
# 或者: from Euclid import adjacency_to_laplacian [as 别名]
def edges_to_laplacian(edges, edge_weights):
    """
    Return a Laplacian matrix given a list of unweighted edges.
    Vertices of the graph are expected to be indexed
    sequentially starting at zero.
    @param edges: a list of connected point index pairs
    @param edge_weights: a list of edge weights
    @return: a numpy array laplacian matrix
    """
    source_indices, sink_indices = zip(*edges)
    all_indices = set(source_indices) | set(sink_indices)
    npoints = max(all_indices) + 1
    if all_indices != set(range(npoints)):
        raise ValueError('the graph is not connected')
    # create the adjacency matrix
    A = np.zeros((npoints, npoints))
    for weight, (i, j) in zip(edge_weights, edges):
        A[i, j] = weight
        A[j, i] = weight
    # create the laplacian matrix
    L = Euclid.adjacency_to_laplacian(A)
    return L

# 需要导入模块: import Euclid [as 别名]
# 或者: from Euclid import adjacency_to_laplacian [as 别名]
def get_response_content(fs):
    M = create_adjacency_matrix(fs.edge_affinity, fs.nvertices)
    if fs.laplacian:
        M = Euclid.adjacency_to_laplacian(M)
    return MatrixUtil.m_to_string(M) + '\n'

# 需要导入模块: import Euclid [as 别名]
# 或者: from Euclid import adjacency_to_laplacian [as 别名]
def get_response_content(fs):
    A = fs.matrix
    L = Euclid.adjacency_to_laplacian(A)
    D = Euclid.laplacian_to_edm(L)
    return MatrixUtil.m_to_string(D) + "\n"