Rコード最適化例:空間点パターンのメディアン
をテンプレートにして作成
[
トップ
] [
新規
|
一覧
|
検索
|
最終更新
|
ヘルプ
]
開始行:
SIZE(20){COLOR(magenta){Rコード最適化例:空間点パターンのメディアン}} [[Rコード最適化のコツと実例集]] に戻る~
[[Q&A]] コーナ投稿記事より
~
空間メジアン(spatial median)と2変量オヤ・メジアン(bivariate Oja median)を計算する。
T(θ) = Σ|xi-θ|とするとき T(θ) を最小にするθが空間メジアン、
T(θ)=ΣA(xi,xj,θ)とするとき T(θ) を最小にするθが2変量オヤ・メジアン
(A(a,b,c)は3点a,b,cを頂点とする三角形の面積)。
## 空間メディアンを optim() で求める
> X <- runif(1000) # テスト用の1000個の点の x 座標
> Y <- runif(1000) # テスト用の1000個の点の y 座標
## 任意の点 T=(T[1], T[2]) から点 (X[i],Y[i]) までの距離和関数
dsum <- function (T) { # 素朴な定義
s <- 0
for (i in 1:length(X)) s <- s + sqrt((T[1]-X[i])^2+(T[2]-Y[i])^2)
s }
> system.time(res <- optim(c(mean(X),mean(Y)), dsum)) # 重心を初期値にとり Nelder-Mead 法で最小化
[1] 0.62 0.00 1.20 0.00 0.00 # 関数評価回数41回
> dsum2 <- function(T) { sum(sqrt((T[1]-X)^2+(T[2]-Y)^2)) } # ベクトル化演算で高速化
> system.time(res2 <- optim(c(mean(X),mean(Y)), dsum2)) # 重心を初期値にとり Nelder-Mead 法で最小化
[1] 0.03 0.00 0.07 0.00 0.00
## 2変量オヤ・メジアンを optim() で求める
## 三点 (X[i],Y[i]),(X[j],Y[j]), (T[1], T[2]) の作る三角形の面積は
## abs((X[i]-T[i])*(Y[j]-T[j]) - (Y[i]-T[i])*(X[j]-T[j]))/2 に注意する
## 三角形面積和を計算する素朴な関数
asum <- function (T) {
asum <- 0
for (i in 1:(length(X)-1)) # 全ての 1 <= i < j <= length(X) に対する和
for (j in (i+1):length(X))
asum <- asum + abs((X[i]-T[1])*(Y[j]-T[2]) - (Y[i]-T[1])*(X[j]-T[2]))
asum/2 } # 毎回割らずに、最後に割っていることを注意
## ベクトルの外積演算を使い高速化
## (X-T[1])%o%(Y-T[2]) は (i,j) 成分が (X[i]-T[1])*(Y[j]-T[2]) である行列
> asum2 <- function (T) {Z <- (X-T[1])%o%(Y-T[2]); sum(abs((Z-t(Z))))/4}
> asum2(c(0,0)) # 検査
[1] 665.7858
> X <- runif(100) # テスト用の100個の点の x 座標
> Y <- runif(100) # テスト用の100個の点の y 座標
> system.time(res <- optim(c(mean(X),mean(Y)), asum)) # 重心を初期値にとり Nelder-Mead 法で最小化
[1] 5.6 0.0 11.1 0.0 0.0 # 関数評価回数 47 回
> system.time(res2 <- optim(c(mean(X),mean(Y)), asum2)) # 重心を初期値にとり Nelder-Mead 法で最小化
[1] 0.09 0.00 0.14 0.00 0.00
> X <- runif(200) # 点の数を200にする
> Y <- runif(200)
> system.time(res2 <- optim(c(mean(X),mean(Y)), asum2)) # 重心を初期値にとり Nelder-Mead 法で最小化
[1] 0.43 0.03 0.86 0.00 0.00 # 関数評価回数47回
> system.time(res <- optim(c(mean(X),mean(Y)), asum)) # 重心を初期値にとり Nelder-Mead 法で最小化
[1] 21.11 0.29 43.64 0.00 0.00
終了行:
SIZE(20){COLOR(magenta){Rコード最適化例:空間点パターンのメディアン}} [[Rコード最適化のコツと実例集]] に戻る~
[[Q&A]] コーナ投稿記事より
~
空間メジアン(spatial median)と2変量オヤ・メジアン(bivariate Oja median)を計算する。
T(θ) = Σ|xi-θ|とするとき T(θ) を最小にするθが空間メジアン、
T(θ)=ΣA(xi,xj,θ)とするとき T(θ) を最小にするθが2変量オヤ・メジアン
(A(a,b,c)は3点a,b,cを頂点とする三角形の面積)。
## 空間メディアンを optim() で求める
> X <- runif(1000) # テスト用の1000個の点の x 座標
> Y <- runif(1000) # テスト用の1000個の点の y 座標
## 任意の点 T=(T[1], T[2]) から点 (X[i],Y[i]) までの距離和関数
dsum <- function (T) { # 素朴な定義
s <- 0
for (i in 1:length(X)) s <- s + sqrt((T[1]-X[i])^2+(T[2]-Y[i])^2)
s }
> system.time(res <- optim(c(mean(X),mean(Y)), dsum)) # 重心を初期値にとり Nelder-Mead 法で最小化
[1] 0.62 0.00 1.20 0.00 0.00 # 関数評価回数41回
> dsum2 <- function(T) { sum(sqrt((T[1]-X)^2+(T[2]-Y)^2)) } # ベクトル化演算で高速化
> system.time(res2 <- optim(c(mean(X),mean(Y)), dsum2)) # 重心を初期値にとり Nelder-Mead 法で最小化
[1] 0.03 0.00 0.07 0.00 0.00
## 2変量オヤ・メジアンを optim() で求める
## 三点 (X[i],Y[i]),(X[j],Y[j]), (T[1], T[2]) の作る三角形の面積は
## abs((X[i]-T[i])*(Y[j]-T[j]) - (Y[i]-T[i])*(X[j]-T[j]))/2 に注意する
## 三角形面積和を計算する素朴な関数
asum <- function (T) {
asum <- 0
for (i in 1:(length(X)-1)) # 全ての 1 <= i < j <= length(X) に対する和
for (j in (i+1):length(X))
asum <- asum + abs((X[i]-T[1])*(Y[j]-T[2]) - (Y[i]-T[1])*(X[j]-T[2]))
asum/2 } # 毎回割らずに、最後に割っていることを注意
## ベクトルの外積演算を使い高速化
## (X-T[1])%o%(Y-T[2]) は (i,j) 成分が (X[i]-T[1])*(Y[j]-T[2]) である行列
> asum2 <- function (T) {Z <- (X-T[1])%o%(Y-T[2]); sum(abs((Z-t(Z))))/4}
> asum2(c(0,0)) # 検査
[1] 665.7858
> X <- runif(100) # テスト用の100個の点の x 座標
> Y <- runif(100) # テスト用の100個の点の y 座標
> system.time(res <- optim(c(mean(X),mean(Y)), asum)) # 重心を初期値にとり Nelder-Mead 法で最小化
[1] 5.6 0.0 11.1 0.0 0.0 # 関数評価回数 47 回
> system.time(res2 <- optim(c(mean(X),mean(Y)), asum2)) # 重心を初期値にとり Nelder-Mead 法で最小化
[1] 0.09 0.00 0.14 0.00 0.00
> X <- runif(200) # 点の数を200にする
> Y <- runif(200)
> system.time(res2 <- optim(c(mean(X),mean(Y)), asum2)) # 重心を初期値にとり Nelder-Mead 法で最小化
[1] 0.43 0.03 0.86 0.00 0.00 # 関数評価回数47回
> system.time(res <- optim(c(mean(X),mean(Y)), asum)) # 重心を初期値にとり Nelder-Mead 法で最小化
[1] 21.11 0.29 43.64 0.00 0.00
ページ名: