미국 전역의 기상 관측소 네트워크에 대한 데이터가 있습니다. 이것은 날짜, 위도, 경도 및 일부 측정 값을 포함하는 데이터 프레임을 제공합니다. 데이터가 하루에 한 번 수집되고 지역 규모의 날씨에 의해 구동된다고 가정합니다 (아니, 우리는 그 논의에 참여하지 않을 것입니다).
동시에 측정 된 값이 시간과 공간에서 어떻게 상관되는지 그래픽으로 보여주고 싶습니다. 저의 목표는 조사중인 가치의 지역적 균질성 (또는 부족)을 보여주는 것입니다.
데이터 세트
우선 매사추세츠 주와 메인 주에있는 여러 역을 가졌습니다. NOAA의 FTP 사이트에서 사용할 수있는 색인 파일에서 위도와 경도로 사이트를 선택했습니다.
바로 하나의 문제가 있습니다. 비슷한 식별자를 가지고 있거나 매우 가까운 사이트가 많이 있습니다. FWIW, 나는 USAF와 WBAN 코드를 모두 사용하여 식별합니다. 메타 데이터를 자세히 살펴보면 좌표와 고도가 다르고 데이터가 한 사이트에서 중지 된 다음 다른 사이트에서 시작한다는 것을 알았습니다. 그래서 나는 더 잘 알지 못하기 때문에 별도의 스테이션으로 취급해야합니다. 이것은 데이터가 서로 매우 가까운 스테이션 쌍을 포함한다는 것을 의미합니다.
예비 분석
나는 월별로 데이터를 그룹화 한 다음 서로 다른 데이터 쌍 사이의 최소 최소 제곱 회귀를 계산하려고했습니다. 그런 다음 스테이션을 연결하는 선으로 모든 쌍 간의 상관 관계를 그립니다 (아래). 선 색상은 OLS 피팅의 R2 값을 보여줍니다. 그런 다음 그림은 1 월, 2 월 등의 30 개 이상의 데이터 포인트가 관심 영역의 다른 스테이션간에 어떻게 상호 연관되는지 보여줍니다.
6 시간마다 데이터 포인트가있는 경우에만 일일 평균이 계산되도록 기본 코드를 작성 했으므로 사이트간에 데이터를 비교할 수 있습니다.
문제
불행히도, 한 플롯에서 이해하기에는 너무 많은 데이터가 있습니다. 선의 크기를 줄이면 해결할 수 없습니다.
네트워크가 너무 복잡해 보이므로 복잡성을 줄이거 나 공간 커널을 적용 할 수있는 방법을 찾아야한다고 생각합니다.
또한 상관 관계를 나타내는 데 가장 적합한 메트릭이 무엇인지 확실하지 않지만 의도 된 (비 기술적 인) 대상의 경우 OLS의 상관 계수가 가장 간단하게 설명 할 수 있습니다. 그래디언트 또는 표준 오류와 같은 다른 정보도 제시해야 할 수도 있습니다.
질문
나는이 분야와 R로 동시에 길을 배우고 있으며 다음에 대한 제안을 부탁드립니다.
- 내가하려는 일의 더 공식적인 이름은 무엇입니까? 더 많은 문헌을 찾을 수있는 유용한 용어가 있습니까? 내 검색에서 공통 응용 프로그램이 무엇인지 빈칸으로 표시됩니다.
- 공간에서 분리 된 여러 데이터 세트 간의 상관 관계를 표시하는 데 더 적합한 방법이 있습니까?
- … 특히 시각적으로 결과를 쉽게 표시 할 수있는 방법은 무엇입니까?
- 이것들은 R로 구현 되었습니까?
- 이러한 접근 방식 중 어느 것이 자동화에 적합합니까?
답변
이 유형의 데이터를 표시하는 몇 가지 옵션이 있다고 생각합니다.
첫 번째 옵션은 “실증적 직교 함수 분석”(EOF) (기후가 아닌 원에서는 “주성분 분석”(PCA)이라고도 함)을 수행하는 것입니다. 귀하의 경우, 이것은 데이터 위치의 상관 매트릭스에서 수행되어야합니다. 예를 들어, 데이터 매트릭스 dat
는 열 차원의 공간 위치이고 행의 측정 된 매개 변수입니다. 따라서 데이터 매트릭스는 각 위치에 대한 시계열로 구성됩니다. 이 prcomp()
기능을 사용하면이 필드와 관련된 주요 구성 요소 또는 지배적 상관 모드를 얻을 수 있습니다.
res <- prcomp(dat, retx = TRUE, center = TRUE, scale = TRUE) # center and scale should be "TRUE" for an analysis of dominant correlation modes)
#res$x and res$rotation will contain the PC modes in the temporal and spatial dimension, respectively.
두 번째 옵션은 관심있는 개별 위치와 관련된 상관 관계를 나타내는 맵을 만드는 것입니다.
C <- cor(dat)
#C[,n] would be the correlation values between the nth location (e.g. dat[,n]) and all other locations.
편집 : 추가 예
다음 예제에서는 갭피 데이터를 사용하지 않지만 DINEOF를 사용하여 보간 한 후 동일한 분석을 데이터 필드에 적용 할 수 있습니다 ( http://menugget.blogspot.de/2012/10/dineof-data-interpolating-empirical.html ) . 아래 예는 다음 데이터 세트 ( http://www.esrl.noaa.gov/psd/gcos_wgsp/Gridded/data.hadslp2.html ) 에서 월별 이상 해수면 압력 데이터의 하위 집합을 사용합니다 .
library(sinkr) # https://github.com/marchtaylor/sinkr
# load data
data(slp)
grd <- slp$grid
time <- slp$date
field <- slp$field
# make anomaly dataset
slp.anom <- fieldAnomaly(field, time)
# EOF/PCA of SLP anom
P <- prcomp(slp.anom, center = TRUE, scale. = TRUE)
expl.var <- P$sdev^2 / sum(P$sdev^2) # explained variance
cum.expl.var <- cumsum(expl.var) # cumulative explained variance
plot(cum.expl.var)
주요 EOF 모드 매핑
# make interpolation
require(akima)
require(maps)
eof.num <- 1
F1 <- interp(x=grd$lon, y=grd$lat, z=P$rotation[,eof.num]) # interpolated spatial EOF mode
png(paste0("EOF_mode", eof.num, ".png"), width=7, height=6, units="in", res=400)
op <- par(ps=10) #settings before layout
layout(matrix(c(1,2), nrow=2, ncol=1, byrow=TRUE), heights=c(4,2), widths=7)
#layout.show(2) # run to see layout; comment out to prevent plotting during .pdf
par(cex=1) # layout has the tendency change par()$cex, so this step is important for control
par(mar=c(4,4,1,1)) # I usually set my margins before each plot
pal <- jetPal
image(F1, col=pal(100))
map("world", add=TRUE, lwd=2)
contour(F1, add=TRUE, col="white")
box()
par(mar=c(4,4,1,1)) # I usually set my margins before each plot
plot(time, P$x[,eof.num], t="l", lwd=1, ylab="", xlab="")
plotRegionCol()
abline(h=0, lwd=2, col=8)
abline(h=seq(par()$yaxp[1], par()$yaxp[2], len=par()$yaxp[3]+1), col="white", lty=3)
abline(v=seq.Date(as.Date("1800-01-01"), as.Date("2100-01-01"), by="10 years"), col="white", lty=3)
box()
lines(time, P$x[,eof.num])
mtext(paste0("EOF ", eof.num, " [expl.var = ", round(expl.var[eof.num]*100), "%]"), side=3, line=1)
par(op)
dev.off() # closes device
상관 관계 맵 작성
loc <- c(-90, 0)
target <- which(grd$lon==loc[1] & grd$lat==loc[2])
COR <- cor(slp.anom)
F1 <- interp(x=grd$lon, y=grd$lat, z=COR[,target]) # interpolated spatial EOF mode
png(paste0("Correlation_map", "_lon", loc[1], "_lat", loc[2], ".png"), width=7, height=5, units="in", res=400)
op <- par(ps=10) #settings before layout
layout(matrix(c(1,2), nrow=2, ncol=1, byrow=TRUE), heights=c(4,1), widths=7)
#layout.show(2) # run to see layout; comment out to prevent plotting during .pdf
par(cex=1) # layout has the tendency change par()$cex, so this step is important for control
par(mar=c(4,4,1,1)) # I usually set my margins before each plot
pal <- colorRampPalette(c("blue", "cyan", "yellow", "red", "yellow", "cyan", "blue"))
ncolors <- 100
breaks <- seq(-1,1,,ncolors+1)
image(F1, col=pal(ncolors), breaks=breaks)
map("world", add=TRUE, lwd=2)
contour(F1, add=TRUE, col="white")
box()
par(mar=c(4,4,0,1)) # I usually set my margins before each plot
imageScale(F1, col=pal(ncolors), breaks=breaks, axis.pos = 1)
mtext("Correlation [R]", side=1, line=2.5)
box()
par(op)
dev.off() # closes device
답변
나는 선 뒤에 명확하게 보이지 않지만 너무 많은 데이터 포인트가있는 것 같습니다.
정확한 스테이션이 아니라 지역적 균질성을 보여주기를 원하기 때문에 먼저 공간적으로 그룹화하는 것이 좋습니다. 예를 들어 “피쉬 넷”으로 오버레이하고 모든 셀에서 (매 순간마다) 평균 측정 값을 계산합니다. 이 평균값을 셀 중심에이 방법으로 배치하면 데이터를 래스터화할 수 있습니다 (또는 선을 겹치지 않으려면 모든 셀의 위도 및 경도도 계산할 수 있음). 또는 관리 단위 내부의 평균을 계산하십시오. 그런 다음이 새로운 평균 “스테이션”의 경우 상관 관계를 계산하고 더 적은 수의 선으로 맵을 그릴 수 있습니다.
이것은 또한 모든 영역을 통과하는 임의의 단일 상관 상관 선을 제거 할 수 있습니다.