原点を頂点とする多次元錐の内側かどうか

• df次元空間を考える
• ある点Aがある
• いくつかの点の集合Xsがある
• 原点とXsの各点を結んだ線によって「籠」のようになった範囲は錐
• Aの方向がこの錐の内部かどうかの判定をしよう
• df個の線形独立なベクトルがあるときそれらを用いて、Aの座標を表現できる
  • そのような表現のときに、すべての係数が0以上であるとき、Aはこの錐の内側(周も含めて)
• Xsにdf個より多い点があるときは、df個ずつの組み合わせにして、そのいずれかが作る錐の内側なら、錐全体の内側である
• Xsにdf個より少ない点があるときは、その線形結合であって、すべての係数が正のときに内側
• 今、「すべての係数が正」で話を進めてきたが、「すべての係数が負」の場合も「逆向きの錐」とすることも可能なので、以下のソースでは、その設定を"bothside={TRUE,FALSE}"でオプションにしてある

df<-3 # 次元
N<-6 # 点の数
library(MCMCpack)
Xs<-matrix(rdirichlet(N,rep(1,df)),nrow=N)
#Xs<-matrix(c(1,0,0,1),nrow=2,byrow=TRUE)

if(N>=df){
#点の数が次元以上ならば、内側になるかどうかを試すべくAを作る
 A<-rnorm(df)
}else{
#点の数が次元未満なら、あえて、線形結合になるようなAを作ってやる(すべてが正なら内側)
 A<-apply(Xs *rnorm(N),2,sum)
}


JudgeInside<-function(A,Xs,bothside=TRUE){# bothside=TRUE は両側の錐の内側判定

 B<-NULL
 df<-length(A) # 次元
 N<-length(Xs[,1]) # 点の数
 inside<-FALSE
 coefs<-NULL
 if(N>=df){
  library(gtools)
  cmb<-combinations(N,df,repeats=FALSE)

  loop<-TRUE
  counter<-1
  vs<-NULL

  while(counter<=length(cmb[,1])&loop){
#選んだ点の組み合わせを座標としてAの位置を表す
   tmp<-solve(t(Xs[cmb[counter,],]),A)
   if(length(tmp[tmp<0])==0 || (bothside & length(tmp[tmp>0])==0) ){# すべての要素が0以上の条件
    inside<-TRUE # 錐の内側
    vs<-Xs[cmb[counter,],]
    coefs<-tmp
    loop<-FALSE
   }
   counter<-counter+1
  }
  if(inside){
   B<-t(vs) %*% coefs
  }

 }else{
  tmp<-solve(t(Xs[,1:N]),A[1:N])

  if(sum((A-B)^2)<10^(-16)& (length(tmp[tmp<0])==0|| (bothside & length(tmp[tmp>0])==0))){
   inside=TRUE
   B<-t(Xs) %*% tmp
   coefs<-tmp
  }

 }
 list(inside=inside,A=A,Xs=Xs,coefs=coefs,B=B)
}

JudgeInside(A,Xs)