2017-01-11 7 views
2

私は、Vincent Zoonekyndによって書かれた2つの関数の助けを借りて次のグラフを作成しました(それらはhereで見つけることができます) )。R:2次元点群の2点間の最短測地線を見つける

Points connected to their 3 neighbouring points

その近傍グラフとそのパラメータ「K」はIsometric Feature Mapping用途、これが何であるかを説明できるようにするために。 "k"は各ポイントが直接接続されているポイント数を指定します。それらの距離はお互いの真理値の距離に過ぎません。任意の点とその(k + 1) - 最も近い点(または遠方の任意の点)との間の距離は「測地線」と呼ばれ、そこに到達するために必要なすべての辺の長さの最小合計です。これは、時には真理距離よりもはるかに長い。これは私の図の点AとBの場合です。

ここで点Aから点Bまでの測地線距離を示す黒線を追加したいと思います。segments()というコマンドについて知っていますが、これはおそらく線を追加するのに最適です。最短経路(測地線距離)はDijkstra's Algorithmであり、パッケージigraphに実装されています。しかし、igraphは自分のグラフを解釈することも、渡す必要のあるポイント(頂点)(およびそれらの座標)を自分で見つけることもできません。ところで、k = 18の場合、すなわち、すべての点が直近の18個の点に直接接続されている場合、AとBの間の測地線距離はちょうど真理値の距離になります。


isomap.incidence.matrix <- function (d, eps=NA, k=NA) { 
    stopifnot(xor(is.na(eps), is.na(k))) 
    d <- as.matrix(d) 
    if(!is.na(eps)) { 
    im <- d <= eps 
    } else { 
    im <- apply(d,1,rank) <= k+1 
    diag(im) <- FALSE 
    } 
    im | t(im) 
} 

plot.graph <- function (im,x,y=NULL, ...) { 
    if(is.null(y)) { 
    y <- x[,2] 
    x <- x[,1] 
    } 
    plot(x,y, ...) 
    k <- which( as.vector(im) ) 
    i <- as.vector(col(im))[ k ] 
    j <- as.vector(row(im))[ k ] 
    segments(x[i], y[i], x[j], y[j], col = "grey") 
} 

z <- seq(1.1,3.7,length=140)*pi 

set.seed(4) 
zz <- rnorm(1:length(z))+z*sin(z) 
zz <- cbind(zz,z*cos(z)*seq(3,1,length=length(z))) 

dist.grafik <- dist(zz) 

pca.grafik <- princomp(zz) 

x11(8, 8) 
par(mar=c(0,0,0,0)) 
plot.graph(isomap.incidence.matrix(dist.grafik, k=3), pca.grafik$scores[,1], pca.grafik$scores[,2], 
      xaxt = "n", yaxt = "n", xlab = "", ylab = "", cex = 1.3) 
legend("topright", inset = 0.02, legend = "k = 3", col = "grey", lty = 1, cex = 1.3) 
segments(x0 = -8.57, y0 = -1.11, x1 = -10.83, y1 = -5.6, col = "black", lwd = 2, lty = "dashed") 
text(x = -8.2, y = -1.4, labels = "A", font = 2, cex = 1.2) 
text(x = -11, y = -5.1, labels = "B", font = 2, cex = 1.2) 
+0

グラフに黒い線を表示する方法がわからないという意味で、またはigraphなしでDijkstraのアルゴリズムをどのように再コード化するかを尋ねるという意味で、ネットワーク関連の問題ですか?それとも、グラフをigraphで解釈させるのかという質問ですか? – probaPerception

+0

私は質問を編集してより明確にしました。 – mattu

答えて

2

次のコードでは、それはあなたのケースである重量とIGRAPHオブジェクト、ノード間のユークリッド距離を作成するためにあなたのデータを使用して、あなたを助けるかもしれません。 次に、加重最短パスが見つけられ、これはsp$vpath[[1]]によって返されます。次の例では、それは私がこのスクリプトを実行するにはmattu

isomap.incidence.matrix <- function (d, eps=NA, k=NA) { 
    stopifnot(xor(is.na(eps), is.na(k))) 
    d <- as.matrix(d) 
    if(!is.na(eps)) { 
    im <- d <= eps 
    } else { 
    im <- apply(d,1,rank) <= k+1 
    diag(im) <- FALSE 
    } 
    im | t(im) 
} 

plot.graph <- function (im,x,y=NULL, ...) { 
    if(is.null(y)) { 
    y <- x[,2] 
    x <- x[,1] 
    } 
    plot(x,y, ...) 
    k <- which( as.vector(im) ) 
    i <- as.vector(col(im))[ k ] 
    j <- as.vector(row(im))[ k ] 
    segments(x[i], y[i], x[j], y[j], col = "grey") 
} 

z <- seq(1.1,3.7,length=100)*pi 

set.seed(4) 
zz <- rnorm(1:length(z))+z*sin(z) 
zz <- cbind(zz,z*cos(z)*seq(3,1,length=length(z))) 

dist.grafik <- as.matrix(dist(zz)) 
pca.grafik <- princomp(zz) 

isomap.resul <- function (d, eps=NA, k=NA) { 
    a <- isomap.incidence.matrix(d, eps, k) 
    b <- dist.grafik 
    res <- a * b 
    return(res) 
} 

a <- graph_from_adjacency_matrix(isomap.resul(dist.grafik, k=3), 
           mode = c("undirected"), weight = TRUE) 
sp <- shortest_paths(a, 5, to = 66, mode = c("out", "all", "in"), 
        weights = NULL, output = c("vpath", "epath", "both"), 
        predecessors = FALSE, inbound.edges = FALSE) 

path <- sp$vpath[[1]] 

x11(8, 8) 
par(mar=c(0,0,0,0)) 
plot.graph(isomap.incidence.matrix(dist.grafik, k=3), pca.grafik$scores[,1], pca.grafik$scores[,2], 
      xaxt = "n", yaxt = "n", xlab = "", ylab = "", cex = 1.3) 
legend("topright", inset = 0.02, legend = "k = 3", col = "grey", lty = 1, cex = 1.3) 
segments(x0 = -8.57, y0 = -1.11, x1 = -10.83, y1 = -5.6, col = "black", lwd = 2, lty = "dashed") 
text(x = -8.2, y = -1.4, labels = "A", font = 2, cex = 1.2) 
text(x = -11, y = -5.1, labels = "B", font = 2, cex = 1.2) 

for(i in 2:length(path)){ 
    aa <- pca.grafik$scores[path[i-1], 1] 
    bb <- pca.grafik$scores[path[i-1], 2] 
    cc <- pca.grafik$scores[path[i], 1] 
    dd <- pca.grafik$scores[path[i], 2] 
    segments(aa, bb, cc , dd, lwd = 2) 
} 

からプロットするために溶液でコードを編集したノード番号5と66 間の最短経路である、あなたは明らかにパッケージigraphを必要としています。

私にとっては、測地距離に応じて最短経路に見えます。 enter image description here

希望します。

+0

さて、ノード番号noから始まります。ノード番号5で終了する。 91はAからBへのパスを与えるでしょう。私は 'sp $ vpath [[1]]' 'path'という名前をつけて、グラフへのパスをもっと読みやすくするためのコードを作りました:' for(i in 2:length(path (パス[i-1]、1)、pca.grafik $スコア[パス[i-1]、2]、pca.grafik $スコア[パス[i]、1] 、pca.grafik $スコア[パス[i]、2]、lwd = 2)} '。しかし、そのパスは最短ではないようです。したがって私は間違いを犯したのか、igraphが最後に行ったのか疑問に思っています。あなたの印象を加えてください。 (私はもちろん、コメントに画像をアップロードすることはできません) – mattu

+0

パスを手作業で少し修正し、この結果を得ました: 'path [1] 5 9 10 13 14 16 18 20 22 24 25 26 29 31 34 36 39 40 44 45 50 47 51 52 53 56 57 59 60 62 64 68 70 72 73 74 76 78 79 80 82 83 85 87 88 90 91 '。これは大丈夫ですが、複製可能なソリューションはありません。 – mattu

+1

私はあなたのソリューションからセグメントをプロットするコードを調整し、私はプロットの画像を追加します。私にとっては、最短経路(測地線上のユークリッド距離という意味で)です。あなたはそうは思わない? – probaPerception

関連する問題