Julien Claude の本のおかげでMorphometrics with R
、その matlab 関数と同じことを行う便利なコードがいくつかあります。
彼は完全な Procrustes 距離を計算する関数をいくつか提供しています。これは、matlab 関数が定義されているように、「重ね合わせられた構成の相同な座標間の距離の 2 乗の合計の平方根 (以前は単位サイズにスケーリングされていた)」として定義されています。
# first, scale the coordinates to unit centroid size, and return both the scaled coords and the centroid size
centsiz<-function(M)
{p<-dim(M)[1]
size<-sqrt(sum(apply(M, 2,var))*(p-1))
list("centroid_size" = size,"scaled" = M/size)}
# second, translate the coords so that its centroid is set at the origin
trans1<-function(M){scale(M,scale=F)}
# third, prepare the fPsup function to perform the full Procrustes superimposition of M1 onto M2. In the output, DF is the Full Procrustes distance between M1 and M2.
fPsup<-function(M1, M2) {
k<-ncol(M1)
Z1<-trans1(centsiz(M1)[[2]])
Z2<-trans1(centsiz(M2)[[2]])
sv<-svd(t(Z2)%*%Z1)
U<-sv$v; V<-sv$u; Delt<-sv$d
sig<-sign(det(t(Z2)%*%Z1))
Delt[k]<-sig*abs(Delt[k]) ; V[,k]<-sig * V[,k]
Gam<-U%*%t(V)
beta<-sum(Delt)
list(Mp1=beta*Z1%*%Gam,Mp2=Z2,rotation=Gam,scale=beta,
DF=sqrt(1-beta^2))}
# test it out...
library(shapes) # so we can use the built-in data
data(gorf.dat) # Female gorilla skull data, 8 landmarks in 2 dimensions, 30 individuals
# calculate procrustes distance for individuals 1 and 2
fPsup(gorf.dat[,,1], gorf.dat[,,2])$DF
[1] 0.0643504
# Claude provides a check with a function that calculates the interlandmark distances between two configurations, which we can then sqrt the sum of to get the matlab-defined procrustes distance.
ild2<-function(M1, M2){sqrt(apply((M1-M2)^2, 1, sum))}
# test it out...
test<-fPsup(gorf.dat[,,1], gorf.dat[,,2])
test$DF
[1] 0.0643504
sqrt(sum(ild2(test$Mp1, test$Mp2)^2))
[1] 0.0643504 # the same
パッケージに固執したいだけならshapes
、リーマン形状距離関数はほぼ同じ結果を計算します。
library(shapes)
riemdist(gorf.dat[,,1], gorf.dat[,,2])
[1] 0.0643949
更新shapes
パッケージの作成者である Ian Dryden と連絡を取りました。彼は、Procrustes の完全な距離を取得するには、 を使用するだけでよいと書いていますsin(riemdist)
。したがって、最初の 2 匹のメスのゴリラ間の完全なプロクラステス距離は次のようになります。
sin(riemdist(gorf.dat[,,1],gorf.dat[,,2]))
[1] 0.0643504
そして、同じことを行う独自の関数を作成したい場合fpdist
:
fpdist<-function(x, y, reflect = FALSE){
sin(riemdist(x,y,reflect=reflect))
}
fpdist(gorf.dat[,,1],gorf.dat[,,2])
[1] 0.0643504
上記で使用されているゴリラ データは 2D ですが、3D データも正常に機能することに注意してください。
library(shapes) # so we can use the built-in data
data(macm.dat) # Male macaque skull data. 7 landmarks in 3 dimensions, 9 individuals
# calculate procrustes distance for macaque individuals 1 and 2
# Claude's method 1
fPsup(macm.dat[,,1], macm.dat[,,2])$DF
[1] 0.1215633
# Claude's method 2
test<-fPsup(macm.dat[,,1], macm.dat[,,2])
sqrt(sum(ild2(test$Mp1, test$Mp2)^2))
[1] 0.1215633
# using the shapes package
fpdist(macm.dat[,,1], macm.dat[,,2])
[1] 0.1215633
それはあなたが求めていたものですか?