Как изменить данные разрешения (или регрессии) в R

Question

У меня есть набор данных, состоящий из lon, lat и среднемесячной переменной (например, температуры или осадков), охватывающий период с 1961 по 1970 год. Набор данных имеет разрешение 0,5 на 0,5 степень долгота/широта и охватывает весь земной шар, и была загружена как файл .NC, который я распакованный данные в R с помощью:

library(ncdf) 
f <- open.ncdf("D:/CRU/cru_ts3.21.1961.1970.tmp.dat.nc") 
A <- get.var.ncdf(nc=f,varid="tmp") 
B <- get.var.ncdf(nc=f,varid="lon") 
C <- get.var.ncdf(nc=f,varid="lat") 
D <- cbind(expand.grid(B, C)) 
E <- expand.grid(A) 

Вспученных сетки (Е) представляет собой таблицу данных, состоящая из 31,104,000 строк переменная и расширенная сетка (D) представляют собой таблицу данных, состоящую из 259 200 строк lon/lat. Если вы умножаете 259 200 * 10 лет * 12 месяцев, вы получите 31,104,000. Следовательно, таблица Е может быть рубят в месячные значения с помощью:

Month <- 1 
Start <- (Month-1)*(259200)+1 
Finish <- (Month*259200) 
G <- E[Start:Finish,] 
H <- expand.grid(G) 
I <- cbind(D,H) 

Поэтому я теперь таблица данных первого месяца (т.е. январь 1961), состоящий из долготы, Lat и переменный. Пример данных приведен ниже:

 lon lat tmp 
49184 -68.25 -55.75 7.5 
49185 -67.75 -55.75 7.6 
49186 -67.25 -55.75 7.6 
49899 -70.75 -55.25 6.8 
49900 -70.25 -55.25 7.0 
49901 -69.75 -55.25 6.9 
49902 -69.25 -55.25 7.1 
49903 -68.75 -55.25 6.8 
49904 -68.25 -55.25 7.6 
49905 -67.75 -55.25 8.2 

Теперь на мой вопрос. Текущее разрешение сетки составляет 0,5 * 0,5 градуса, и я хотел бы «восстановить» данные, чтобы разрешение составляло 0,25 * 0,25 градуса. Я не хочу делать что-либо особенно умное с данными, поэтому я просто хочу, чтобы сетка 0,25 использовала значение сетки 0,5, в которой она сидит, т.е. каждая сетка 0,5 * 0,5 содержит 4 0,25 * 0,25 сетки, и я просто хочу, чтобы 4 0,25 * 0,25 сетки имеют то же значение, что и сетка 0,5 * 0,5.

Я посмотрел растровый, но, похоже, ничего не могу с этим поделать.

Answer 1

Вот как это сделать, используя plyr::ddply(). Вероятно, это будет немного медленнее для вашего размера стола, в зависимости от того, как часто вы хотите переустановить сетку. У меня будет думать о пути, чтобы сделать это с data.table, который должен быть быстрее:

require(plyr) 
# make your data frame 
I<-data.frame(lat=seq(0.5,1000,0.5),lon=1,tmp=sample(1:100,2000,replace=T)) 

# make an adjustment grid 
k<-expand.grid(c(0,0.25),c(0,0.25),0) 

# use plyr:ddply() to expand out each entry into the correponding 4 entries 
new_I<-ddply(I,.(lat,lon),function(x)as.list(x)+k) 
colnames(new_I)<-c("lat","lon","newlat","newlon","tmp") 

head(new_I) 

    lat lon newlat newlon tmp 
1 0.5 1 0.50 1.00 64 
2 0.5 1 0.75 1.00 64 
3 0.5 1 0.50 1.25 64 
4 0.5 1 0.75 1.25 64 
5 1.0 1 1.00 1.00 31 
6 1.0 1 1.25 1.00 31 

На самом деле думать об этом, вот лучший способ с точки зрения времени (хотя это немного рубить , и дает вам меньше контроля за дополнительной обработкой данных, которую вы, возможно, захотите сделать в будущем), но он занимает 6.5 секунд для строк 2m >> 8M.

# make your data frame 
I<-data.frame(lat=seq(0.5,1000000,0.5),lon=1,tmp=sample(1:100,2000000,replace=T)) 

# make an adjustment vector 
v<-rep(0.25,times=2000000) 

# make 3 new tables, apply the vector appropriately, and rbind 
I_latshift<-I 
I_lonshift<-I 
I_bothshift<-I 

I_latshift$lat<-I_latshift$lat+v 
I_lonshift$lon<-I_lonshift$lon+v 
I_bothshift$lat<-I_bothshift$lat+v 
I_bothshift$lon<-I_bothshift$lon+v 

I<-rbind(I,I_bothshift,I_latshift,I_lonshift) 

# sort it for neatness 
I<-I[with(I, order(lat, lon)), ] 


head(I) 

     lat lon tmp 
1  0.50 1.00 3 
6000001 0.50 1.25 3 
4000001 0.75 1.00 3 
2000001 0.75 1.25 3 
2  1.00 1.00 88 
6000002 1.00 1.25 88 

Answer 2

Существует решение в пакете R raster. Он идет следующим образом:

library("ncdf4") 
library("raster") 
nc <- nc_open("my_file.nc") 
lon <- ncvar_get(nc, "lon") 
lat <- ncvar_get(nc, "lat") 
time <- ncvar_get(nc, "time") 
dname <- "pre"  ## pre for the short name of precpitation 
nlon <- dim(lon) 
nlat <- dim(lat) 
nt <- dim(time) 
lonlat <- expand.grid(lon, lat) # make grid of given longitude and latitude 
pr.array <- ncvar_get(nc, dname) 
dlname <- ncatt_get(nc, dname, "long_name") 
dunits <- ncatt_get(nc, dname, "units") 
fillvalue <- ncatt_get(nc, dname, "_FillValue") 

pr.vec.long <- as.vector(pr.array) 
pr.mat <- matrix(pr.vec.long, nrow = nlon * nlat, ncol = nt) 
pr.df <- data.frame(cbind(lonlat, pr.mat)) 

pr_c <- pr.df[ ,-c(1:2)] 
### Specific region have been clipped out from global datafile by 
## selecting lon and lat range and extract regridded data at 1lon 1lat 
## resolution. 

x0 <- seq(67.5, 98.5, by = 1) ## choose different resolution, eg. by = 0.5 
y0 <- seq(6.5, 37.5, by = 1) 


m <- cbind(x0, y0) 
m <- as.data.frame(m) 
s <- rasterFromXYZ(m) 
pts <- expand.grid(x0, y0) 
pos <- pr.df[ ,c(1:2)] 
l_pr <- apply(pr_c, 2, function(x) cbind(pos, x)) 
colnm = c("x","y","z") 
for (j in seq_along(l_pr)){ 
    colnames(l_pr[[j]]) <- colnm 
} 

pr_rstr <- lapply(l_pr, function(x) rasterFromXYZ(x)) 
## Use resample command to regrid the data, here nearest neighbor method can also be chosen by setting method = "ngb" 
pr_bn <- lapply(pr_rstr, function(x) resample(x, s, method = "bilinear")) 
pr_extr <- lapply(pr_bn, function(x) extract(x, pts)) 
df_pr <- do.call("cbind", lapply(pr_extr, data.frame)) 
## write dataframe in csv format 
write.csv(df_pr, "my_data_regridded_1.csv") 

Я надеюсь, что это послужит цели.

Answer 3

Это не решение R, а просто указать, что вы можете легко использовать CDO для восстановления файлов netcdf из командной строки в среде Linux/MAC OS. Из вашего описания, это звучит так, как будто вы хотите использовать ближайший сосед интерполяцию, который для 0.25degree регулярной сетки будет

cdo remapnn,r1440x720 in.nc out.nc 

Однако, вы также можете использовать первый или второй порядок консервативного переназначение. Например, для первого порядка:

cdo remapcon,r1440x720 in.nc out.nc 

Вы можете прочитать в regridded поле в R таким же образом, вы сейчас делаете.