diff --git a/DESCRIPTION b/DESCRIPTION index 86d6601..0a87c9b 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -26,6 +26,6 @@ Suggests: LinkingTo: Rcpp, RcppArmadillo -RoxygenNote: 7.2.3 +RoxygenNote: 7.3.0 Roxygen: list(markdown = TRUE) VignetteBuilder: knitr diff --git a/NAMESPACE b/NAMESPACE index a572513..f3871ec 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -4,6 +4,7 @@ export(annotate_circle) export(draw_circle) export(layout_as_backbone) export(layout_as_dynamic) +export(layout_as_metromap) export(layout_as_multilevel) export(layout_igraph_backbone) export(layout_igraph_centrality) diff --git a/NEWS.md b/NEWS.md index 0536860..e895f58 100644 --- a/NEWS.md +++ b/NEWS.md @@ -3,6 +3,7 @@ * `layout_with_constrained_stress()` and `layout_with_constrained_stress3D()` work for disconnected graphs * internal code refactoring +* added `layout_as_metromap()` # graphlayouts 1.0.2 diff --git a/R/RcppExports.R b/R/RcppExports.R index 9d1b7a0..7f4a637 100644 --- a/R/RcppExports.R +++ b/R/RcppExports.R @@ -17,6 +17,30 @@ constrained_stress_major3D <- function(y, dim, W, D, iter, tol) { .Call(`_graphlayouts_constrained_stress_major3D`, y, dim, W, D, iter, tol) } +criterion_angular_resolution <- function(adj, xy) { + .Call(`_graphlayouts_criterion_angular_resolution`, adj, xy) +} + +criterion_edge_length <- function(el, xy, lg) { + .Call(`_graphlayouts_criterion_edge_length`, el, xy, lg) +} + +criterion_balanced_edge_length <- function(adj_deg2, xy) { + .Call(`_graphlayouts_criterion_balanced_edge_length`, adj_deg2, xy) +} + +criterion_line_straightness <- function() { + .Call(`_graphlayouts_criterion_line_straightness`) +} + +criterion_octilinearity <- function(el, xy) { + .Call(`_graphlayouts_criterion_octilinearity`, el, xy) +} + +layout_as_metro_iter <- function(adj, el, adj_deg2, xy, bbox, l, gr, w, bsize) { + .Call(`_graphlayouts_layout_as_metro_iter`, adj, el, adj_deg2, xy, bbox, l, gr, w, bsize) +} + reweighting <- function(el, N_ranks) { .Call(`_graphlayouts_reweighting`, el, N_ranks) } diff --git a/R/data-examples.R b/R/data-examples.R index ff7feb9..92dad25 100644 --- a/R/data-examples.R +++ b/R/data-examples.R @@ -3,3 +3,11 @@ #' @format igraph object "multilvl_ex" + +#' Subway network of Berlin +#' +#' A dataset containing the subway network of Berlin +#' @format igraph object +#' @references +#' Kujala, Rainer, et al. "A collection of public transport network data sets for 25 cities." Scientific data 5 (2018): 180089. +"metro_berlin" diff --git a/R/metro_multicriteria.R b/R/metro_multicriteria.R new file mode 100644 index 0000000..7d59c84 --- /dev/null +++ b/R/metro_multicriteria.R @@ -0,0 +1,77 @@ +#' @title Metro Map Layout +#' @description Metro map layout based on multicriteria optimization +#' @param object original graph +#' @param xy initial layout of the original graph +#' @param l desired multiple of grid point spacing. (l*gr determines desired edge length) +#' @param gr grid spacing. (l*gr determines desired edge length) +#' @param w weight vector for criteria (see details) +#' @param bsize number of grid points a station can move away rom its original position +#' @details The function optimizes the following five criteria using a hill climbing algorithm: +#' - *Angular Resolution Criterion*: The angles of incident edges at each station should be maximized, because if there is only a small angle between any two adjacent edges, then it can become difficult to distinguish between them +#' - *Edge Length Criterion*: The edge lengths across the whole map should be approximately equal to ensure regular spacing between stations. It is based on the preferred multiple, l, of the grid spacing, g. The purpose of the criterion is to penalize edges that are longer than or shorter than lg. +#' - *Balanced Edge Length Criterion*: The length of edges incident to a particular station should be similar +#' - *Line Straightness Criterion*: (not yet implemented) Edges that form part of a line should, where possible, be co-linear either side of each station that the line passes through +#' - *Octiinearity Criterion*: Each edge should be drawn horizontally, vertically, or diagonally at 45 degree, so we penalize edges that are not at a desired angle +#' @return new coordinates for stations +#' @references +#' Stott, Jonathan, et al. "Automatic metro map layout using multicriteria optimization." IEEE Transactions on Visualization and Computer Graphics 17.1 (2010): 101-114. +#' @author David Schoch +#' @examples +#' # the algorithm has problems with parallel edges +#' library(igraph) +#' g <- simplify(metro_berlin) +#' xy <- cbind(V(g)$lon, V(g)$lat) * 100 +#' +#' # the algorithm is not very stable. try playing with the parameters +#' xy_new <- layout_as_metromap(g, xy, l = 2, gr = 0.5, w = c(100, 100, 1, 1, 100), bsize = 35) +#' @export +layout_as_metromap <- function(object, xy, l = 2, gr = 0.0025, w = rep(1, 5), bsize = 5) { + adj <- as_adj_list1(object) + adj <- lapply(adj, function(x) x - 1) + adj_deg2 <- adj[unlist(lapply(adj, length)) == 2] + el <- igraph::get.edgelist(object, FALSE) - 1 + + xy <- snap_to_grid(xy, gr) + + bbox <- station_bbox(xy, bsize, gr) + + xy_new <- layout_as_metro_iter(adj, el, adj_deg2, xy, bbox, l, gr, w, bsize) + xy_new +} + +# helper ---- +snap_to_grid <- function(xy, gr) { + xmin <- min(xy[, 1]) + xmax <- max(xy[, 1]) + ymin <- min(xy[, 2]) + ymax <- max(xy[, 2]) + + deltax <- seq(xmin - 4 * gr, xmax + 4 * gr, by = gr) + deltay <- seq(ymin - 4 * gr, ymax + 4 * gr, by = gr) + + xdiff <- outer(xy[, 1], deltax, function(x, y) abs(x - y)) + ydiff <- outer(xy[, 2], deltay, function(x, y) abs(x - y)) + + xy_new <- cbind(deltax[apply(xdiff, 1, which.min)], deltay[apply(ydiff, 1, which.min)]) + dups <- duplicated(xy_new) + while (any(dups)) { + xy_new[which(dups), ] <- xy_new[which(dups), ] + c(sample(c(1, -1), 1) * gr, sample(c(1, -1), 1) * gr) + dups <- duplicated(xy_new) + } + xy_new +} + +station_bbox <- function(xy, bsize, gr) { + cbind(xy - bsize * gr, xy + bsize * gr) +} + +as_adj_list1 <- function(g) { + n <- igraph::vcount(g) + lapply(1:n, function(i) { + x <- g[[i]][[1]] + attr(x, "env") <- NULL + attr(x, "graph") <- NULL + class(x) <- NULL + x + }) +} diff --git a/data/metro_berlin.rda b/data/metro_berlin.rda new file mode 100644 index 0000000..8168056 Binary files /dev/null and b/data/metro_berlin.rda differ diff --git a/data/multilvl_ex.rda b/data/multilvl_ex.rda index a95a62b..5e2c864 100644 Binary files a/data/multilvl_ex.rda and b/data/multilvl_ex.rda differ diff --git a/man/layout_as_metromap.Rd b/man/layout_as_metromap.Rd new file mode 100644 index 0000000..5e89e46 --- /dev/null +++ b/man/layout_as_metromap.Rd @@ -0,0 +1,52 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/metro_multicriteria.R +\name{layout_as_metromap} +\alias{layout_as_metromap} +\title{Metro Map Layout} +\usage{ +layout_as_metromap(object, xy, l = 2, gr = 0.0025, w = rep(1, 5), bsize = 5) +} +\arguments{ +\item{object}{original graph} + +\item{xy}{initial layout of the original graph} + +\item{l}{desired multiple of grid point spacing. (l*gr determines desired edge length)} + +\item{gr}{grid spacing. (l*gr determines desired edge length)} + +\item{w}{weight vector for criteria (see details)} + +\item{bsize}{number of grid points a station can move away rom its original position} +} +\value{ +new coordinates for stations +} +\description{ +Metro map layout based on multicriteria optimization +} +\details{ +The function optimizes the following five criteria using a hill climbing algorithm: +\itemize{ +\item \emph{Angular Resolution Criterion}: The angles of incident edges at each station should be maximized, because if there is only a small angle between any two adjacent edges, then it can become difficult to distinguish between them +\item \emph{Edge Length Criterion}: The edge lengths across the whole map should be approximately equal to ensure regular spacing between stations. It is based on the preferred multiple, l, of the grid spacing, g. The purpose of the criterion is to penalize edges that are longer than or shorter than lg. +\item \emph{Balanced Edge Length Criterion}: The length of edges incident to a particular station should be similar +\item \emph{Line Straightness Criterion}: (not yet implemented) Edges that form part of a line should, where possible, be co-linear either side of each station that the line passes through +\item \emph{Octiinearity Criterion}: Each edge should be drawn horizontally, vertically, or diagonally at 45 degree, so we penalize edges that are not at a desired angle +} +} +\examples{ +# the algorithm has problems with parallel edges +library(igraph) +g <- simplify(metro_berlin) +xy <- cbind(V(g)$lon, V(g)$lat) * 100 + +# the algorithm is not very stable. try playing with the parameters +xy_new <- layout_as_metromap(g, xy, l = 2, gr = 0.5, w = c(100, 100, 1, 1, 100), bsize = 35) +} +\references{ +Stott, Jonathan, et al. "Automatic metro map layout using multicriteria optimization." IEEE Transactions on Visualization and Computer Graphics 17.1 (2010): 101-114. +} +\author{ +David Schoch +} diff --git a/man/metro_berlin.Rd b/man/metro_berlin.Rd new file mode 100644 index 0000000..86e8fde --- /dev/null +++ b/man/metro_berlin.Rd @@ -0,0 +1,19 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/data-examples.R +\docType{data} +\name{metro_berlin} +\alias{metro_berlin} +\title{Subway network of Berlin} +\format{ +igraph object +} +\usage{ +metro_berlin +} +\description{ +A dataset containing the subway network of Berlin +} +\references{ +Kujala, Rainer, et al. "A collection of public transport network data sets for 25 cities." Scientific data 5 (2018): 180089. +} +\keyword{datasets} diff --git a/src/RcppExports.cpp b/src/RcppExports.cpp index f6e912a..8193d3d 100644 --- a/src/RcppExports.cpp +++ b/src/RcppExports.cpp @@ -69,6 +69,84 @@ BEGIN_RCPP return rcpp_result_gen; END_RCPP } +// criterion_angular_resolution +double criterion_angular_resolution(List adj, NumericMatrix xy); +RcppExport SEXP _graphlayouts_criterion_angular_resolution(SEXP adjSEXP, SEXP xySEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< List >::type adj(adjSEXP); + Rcpp::traits::input_parameter< NumericMatrix >::type xy(xySEXP); + rcpp_result_gen = Rcpp::wrap(criterion_angular_resolution(adj, xy)); + return rcpp_result_gen; +END_RCPP +} +// criterion_edge_length +double criterion_edge_length(IntegerMatrix el, NumericMatrix xy, double lg); +RcppExport SEXP _graphlayouts_criterion_edge_length(SEXP elSEXP, SEXP xySEXP, SEXP lgSEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< IntegerMatrix >::type el(elSEXP); + Rcpp::traits::input_parameter< NumericMatrix >::type xy(xySEXP); + Rcpp::traits::input_parameter< double >::type lg(lgSEXP); + rcpp_result_gen = Rcpp::wrap(criterion_edge_length(el, xy, lg)); + return rcpp_result_gen; +END_RCPP +} +// criterion_balanced_edge_length +double criterion_balanced_edge_length(List adj_deg2, NumericMatrix xy); +RcppExport SEXP _graphlayouts_criterion_balanced_edge_length(SEXP adj_deg2SEXP, SEXP xySEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< List >::type adj_deg2(adj_deg2SEXP); + Rcpp::traits::input_parameter< NumericMatrix >::type xy(xySEXP); + rcpp_result_gen = Rcpp::wrap(criterion_balanced_edge_length(adj_deg2, xy)); + return rcpp_result_gen; +END_RCPP +} +// criterion_line_straightness +double criterion_line_straightness(); +RcppExport SEXP _graphlayouts_criterion_line_straightness() { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + rcpp_result_gen = Rcpp::wrap(criterion_line_straightness()); + return rcpp_result_gen; +END_RCPP +} +// criterion_octilinearity +double criterion_octilinearity(IntegerMatrix el, NumericMatrix xy); +RcppExport SEXP _graphlayouts_criterion_octilinearity(SEXP elSEXP, SEXP xySEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< IntegerMatrix >::type el(elSEXP); + Rcpp::traits::input_parameter< NumericMatrix >::type xy(xySEXP); + rcpp_result_gen = Rcpp::wrap(criterion_octilinearity(el, xy)); + return rcpp_result_gen; +END_RCPP +} +// layout_as_metro_iter +NumericMatrix layout_as_metro_iter(List adj, IntegerMatrix el, List adj_deg2, NumericMatrix xy, NumericMatrix bbox, double l, double gr, NumericVector w, double bsize); +RcppExport SEXP _graphlayouts_layout_as_metro_iter(SEXP adjSEXP, SEXP elSEXP, SEXP adj_deg2SEXP, SEXP xySEXP, SEXP bboxSEXP, SEXP lSEXP, SEXP grSEXP, SEXP wSEXP, SEXP bsizeSEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< List >::type adj(adjSEXP); + Rcpp::traits::input_parameter< IntegerMatrix >::type el(elSEXP); + Rcpp::traits::input_parameter< List >::type adj_deg2(adj_deg2SEXP); + Rcpp::traits::input_parameter< NumericMatrix >::type xy(xySEXP); + Rcpp::traits::input_parameter< NumericMatrix >::type bbox(bboxSEXP); + Rcpp::traits::input_parameter< double >::type l(lSEXP); + Rcpp::traits::input_parameter< double >::type gr(grSEXP); + Rcpp::traits::input_parameter< NumericVector >::type w(wSEXP); + Rcpp::traits::input_parameter< double >::type bsize(bsizeSEXP); + rcpp_result_gen = Rcpp::wrap(layout_as_metro_iter(adj, el, adj_deg2, xy, bbox, l, gr, w, bsize)); + return rcpp_result_gen; +END_RCPP +} // reweighting NumericVector reweighting(IntegerMatrix el, List N_ranks); RcppExport SEXP _graphlayouts_reweighting(SEXP elSEXP, SEXP N_ranksSEXP) { @@ -191,6 +269,12 @@ static const R_CallMethodDef CallEntries[] = { {"_graphlayouts_constrained_stress_major", (DL_FUNC) &_graphlayouts_constrained_stress_major, 6}, {"_graphlayouts_constrained_stress3D", (DL_FUNC) &_graphlayouts_constrained_stress3D, 3}, {"_graphlayouts_constrained_stress_major3D", (DL_FUNC) &_graphlayouts_constrained_stress_major3D, 6}, + {"_graphlayouts_criterion_angular_resolution", (DL_FUNC) &_graphlayouts_criterion_angular_resolution, 2}, + {"_graphlayouts_criterion_edge_length", (DL_FUNC) &_graphlayouts_criterion_edge_length, 3}, + {"_graphlayouts_criterion_balanced_edge_length", (DL_FUNC) &_graphlayouts_criterion_balanced_edge_length, 2}, + {"_graphlayouts_criterion_line_straightness", (DL_FUNC) &_graphlayouts_criterion_line_straightness, 0}, + {"_graphlayouts_criterion_octilinearity", (DL_FUNC) &_graphlayouts_criterion_octilinearity, 2}, + {"_graphlayouts_layout_as_metro_iter", (DL_FUNC) &_graphlayouts_layout_as_metro_iter, 9}, {"_graphlayouts_reweighting", (DL_FUNC) &_graphlayouts_reweighting, 2}, {"_graphlayouts_sparseStress", (DL_FUNC) &_graphlayouts_sparseStress, 6}, {"_graphlayouts_stress", (DL_FUNC) &_graphlayouts_stress, 3}, diff --git a/src/metroLayout.cpp b/src/metroLayout.cpp new file mode 100644 index 0000000..6814acf --- /dev/null +++ b/src/metroLayout.cpp @@ -0,0 +1,185 @@ +#include +using namespace Rcpp; + +double angle_between_edges(NumericVector pvec,NumericVector qvec) { + if(pvec[0]==qvec[0] && pvec[1]==qvec[1]){ + return 0.0; + } + double dot_pq = pvec[0] * qvec[0] + pvec[1] * qvec[1]; + double mag_pq = sqrt(pvec[0]*pvec[0]+pvec[1]*pvec[1]) * sqrt(qvec[0]*qvec[0]+qvec[1]*qvec[1]); + + if(dot_pq/mag_pq < -0.99){ + return M_PI; + } + if(dot_pq/mag_pq > 0.99){ + return 0.0; + } + + return acos(dot_pq/mag_pq); +} + +// [[Rcpp::export]] +double criterion_angular_resolution(List adj,NumericMatrix xy){ + int n = adj.length(); + double crit=0; + double elen; + double angle; + + for(int i=0;i=bbox(v,0)) && (x<=bbox(v,2)) && (y>=bbox(v,1)) && (y<=bbox(v,3))){ + xy(v,0) = x; + xy(v,1) = y; + tmp_crit = criterion_sum(adj,el,adj_deg2,xy,lg,w); + if(tmp_crit < cur_crit){ + cur_crit = tmp_crit; + running = true; + xbest = xy(v,0); + ybest = xy(v,1); + } + } + } + xy(v,0) = xbest; + xy(v,1) = ybest; + } + // Rcout << "run: " << k << " min: " << cur_crit << std::endl; + k +=1; + } + return xy; +} diff --git a/src/stress.cpp b/src/stress.cpp index e88c685..83d05bc 100644 --- a/src/stress.cpp +++ b/src/stress.cpp @@ -2,129 +2,119 @@ using namespace Rcpp; // [[Rcpp::export]] -double stress(NumericMatrix x, NumericMatrix W, NumericMatrix D){ - double fct=0; - int n=x.nrow(); - for(int i=0;i<(n-1);++i){ - for(int j=(i+1);j0.00001){ - xnew(i,0) += W(i,j)*(x(j,0)+D(i,j)*(x(i,0)-x(j,0))/denom); - xnew(i,1) += W(i,j)*(x(j,1)+D(i,j)*(x(i,1)-x(j,1))/denom); + NumericMatrix x(clone(y)); + + NumericVector wsum = rowSums(W); + + double stress_old = stress(x, W, D); + NumericMatrix xnew(n, 2); // out or in? + for (int k = 0; k < iter; ++k) { + std::fill(xnew.begin(), xnew.end(), 0); + for (int i = 0; i < n; ++i) { + for (int j = 0; j < n; ++j) { + if (i != j) { + double denom = sqrt((x(i, 0) - x(j, 0)) * (x(i, 0) - x(j, 0)) + + (x(i, 1) - x(j, 1)) * (x(i, 1) - x(j, 1))); + if (denom > 0.00001) { + xnew(i, 0) += + W(i, j) * (x(j, 0) + D(i, j) * (x(i, 0) - x(j, 0)) / denom); + xnew(i, 1) += + W(i, j) * (x(j, 1) + D(i, j) * (x(i, 1) - x(j, 1)) / denom); } } } - xnew(i,0) = xnew(i,0)/wsum[i]; - xnew(i,1) = xnew(i,1)/wsum[i]; + xnew(i, 0) = xnew(i, 0) / wsum[i]; + xnew(i, 1) = xnew(i, 1) / wsum[i]; } - double stress_new = stress(xnew,W,D); - double eps = (stress_old-stress_new)/stress_old; + double stress_new = stress(xnew, W, D); + double eps = (stress_old - stress_new) / stress_old; - if(eps<= tol){ + if (eps <= tol) { break; } - stress_old=stress_new; - for(int i=0;i0.00001){ - denom = 1/denom; - } else{ + for (int k = 0; k < iter; ++k) { + NumericMatrix xnew(n, 2); // out or in? + for (int i = 0; i < n; ++i) { + for (int j = 0; j < n; ++j) { + if (i != j) { + double denom = sqrt((x(i, 0) - x(j, 0)) * (x(i, 0) - x(j, 0)) + + (x(i, 1) - x(j, 1)) * (x(i, 1) - x(j, 1))); + if (denom > 0.00001) { + denom = 1 / denom; + } else { denom = 0; } - xnew(i,0) += ((1-t)*W(i,j)+t*Z(i,j))*(x(j,0)+D(i,j)*(x(i,0)-x(j,0))*denom); - xnew(i,1) += ((1-t)*W(i,j)+t*Z(i,j))*(x(j,1)+D(i,j)*(x(i,1)-x(j,1))*denom); + xnew(i, 0) += ((1 - t) * W(i, j) + t * Z(i, j)) * + (x(j, 0) + D(i, j) * (x(i, 0) - x(j, 0)) * denom); + xnew(i, 1) += ((1 - t) * W(i, j) + t * Z(i, j)) * + (x(j, 1) + D(i, j) * (x(i, 1) - x(j, 1)) * denom); } } - xnew(i,0) = xnew(i,0)/((1-t)*wsum[i] + t*zsum[i]); - xnew(i,1) = xnew(i,1)/((1-t)*wsum[i] + t*zsum[i]); + xnew(i, 0) = xnew(i, 0) / ((1 - t) * wsum[i] + t * zsum[i]); + xnew(i, 1) = xnew(i, 1) / ((1 - t) * wsum[i] + t * zsum[i]); } - double stress_newW = stress(xnew,W,D); - double stress_newZ = stress(xnew,Z,D); - double stress_new = (1-t)*stress_newW + t*stress_newZ; + double stress_newW = stress(xnew, W, D); + double stress_newZ = stress(xnew, Z, D); + double stress_new = (1 - t) * stress_newW + t * stress_newZ; - for(int i=0;i