R 'ks'패키지를 사용하여 중복 kdes 내의 데이터 포인트 추출하기

Question

R에서 ks 패키지를 사용하고 데이터 포인트가 중복 된 2 차원 커널 윤곽선의 영역에 속하는 위치를 결정하고 싶습니다. (두 개의 다른 종의 UD를 비교하고 있습니다. 가정 범위). 아래 예가 있습니다 (http://www.rdocumentation.org/packages/ks/versions/1.5.3/topics/plot.kde?에서 수정되었습니다).

생성하려고하는 것은 fhatx의 윤곽선 내에있는 모든 y 점의 목록입니다 (예 : 검은 선 안의 노란색 점). 그리고 그 반대도 마찬가지입니다. 나는 x 좌표의리스트를 윤곽선 안에 넣고 싶습니다.

library(ks) 
x <- rmvnorm.mixt(n=100, mus=c(0,0), Sigmas=diag(2), props=1) 
Hx <- Hpi(x) 
fhatx <- kde(x=x, H=Hx) 
y <- rmvnorm.mixt(n=100, mus=c(0.5,0.5), Sigmas=0.5*diag(2), props=1) 
Hy <- Hpi(y) 
fhaty <- kde(x=y, H=Hy) 
contourLevels(fhatx, cont=c(75, 50, 25)) 
contourSizes(fhatx, cont=25, approx=TRUE) 
plot(fhatx, cont=c(50,95), drawpoints=TRUE) 
plot(fhaty, cont=c(50,95), col=3, drawpoints=TRUE,col.pt="yellow", add=TRUE) 

Answer 1

kde의 출력은 래스터로 변환 될 수 후 거기에서 rasterToPolygons는 함수를 사용하여 임의의 형상을 추출 할 수있다. 포인트가 sp 패키지로 인식되는 형식으로 변환되면 gIntersection 함수를 사용하여 공간 객체 간의 모든 종류의 교차를 볼 수 있습니다.

두 개의 SpatialPoints 객체 x.inY 및 y.inX은 50 % 윤곽선 내에 포함 된 x 점을 포함하며 그 반대도 마찬가지입니다. 이 점의 좌표는 coordinates(...)을 사용하여 배열에서 추출 할 수 있습니다.

이것은 아마도 가장 우아한 해결책은 아니지만 kde 함수로 출시 된 배열이 너무 크지 않으면 잘 작동합니다.

이 정보가 도움이되기를 바랍니다.

STEP 1 : 래스터 객체 KDE로부터의 출력을 변환

# for the x kde 
arrayX <- expand.grid(list(fhatx$eval.points[[1]],fhatx$eval.points[[2]])) 
arrayX$z <- as.vector(fhatx$estimate) 
rasterX <- rasterFromXYZ(arrayX) 
# for the y kde 
arrayY <- expand.grid(list(fhaty$eval.points[[1]],fhaty$eval.points[[2]])) 
arrayY$z <- as.vector(fhaty$estimate) 
rasterY <- rasterFromXYZ(arrayY) 

STEP 2 : 다음 윤곽 월 물론 1. 50 컨투어 내의 모든 셀을 변환하여 0에서 100 사이의 래스터 크기를 조정 (95) 또는 다른 값으로 변경

#for raster x 
rasterX <- rasterX*100/rasterX@data@max 
rasterX[rasterX[]<=50,] <- NA 
rasterX[rasterX[]>50,] <- 1 
#[enter image description here][1]for raster y 
rasterY <- rasterY*100/rasterY@data@max 
rasterY[rasterY[]<=50,] <- NA 
rasterY[rasterY[]>50,] <- 1 

STEP 3 : 폴리곤 50 % 컨투어

polyX50<-rasterToPolygons(rasterX, n=16, na.rm=T, digits=4, dissolve=T) 
polyY50<-rasterToPolygons(rasterY, n=16, na.rm=T, digits=4, dissolve=T) 

대응 추출 691,363,210

4 단계 : 하나의 다각형 또는 다른

#x points falling in fhatx 50 contour 
x.inY <- gIntersection(x.points, polyY50) 
#y points falling in fhatx 50 contour 
y.inX <- gIntersection(y.points, polyX50) 

플롯

와 교차하는 지점을 찾습니다

x.points <- SpatialPoints(x) 
y.points <- SpatialPoints(y) 

STEP 5 SP 라이브러리를 사용하기 위해 공간 객체에 포인트를 변환

par(mfrow=c(1,2)) 
plot(fhatx, cont=c(50,95), col="red") 
plot(fhaty, cont=c(50,95), col="green",add=TRUE) 
plot(x.points, col="red", add=T) 
plot(y.points, col="green", add=T) 

plot(fhatx, cont=c(50,95), col="red") 
plot(fhaty, cont=c(50,95), col="green",add=TRUE) 
plot(x.inY, col="red", add=T) 
plot(y.inX, col="green", add=T)