Welcome to OStack Knowledge Sharing Community for programmer and developer-Open, Learning and Share
Welcome To Ask or Share your Answers For Others

Categories

0 votes
3.7k views
in Technique[技术] by (71.8m points)

grid - Project rasters in R correctly and retrieve their data

I am desperately trying to use two different rasters for getting data and then run regressions with that. I will first try to present what I want (need) to do and then present the problem(s).

Exercise: I am trying to detect the relationship between historical temperature abnormalities and GDP. I did it already on country level, but since my results were unsatisfying, I am redoing the analysis with gridded data by using night-light data as proxy for GDP. I am absolutely new to spatial data analysis and now I am struggling with (probably) an easy solution.

Data:
-Climate ("wilmott" below): https://climatedataguide.ucar.edu/climate-data/global-land-precipitation-and-temperature-willmott-matsuura-university-delaware
-Night-lights:https://ngdc.noaa.gov/eog/dmsp/downloadV4composites.html

Code: Loading both data sources and dividing the climate data brick in individual layers (below I used January 1901 as example, summarizing the layers over the years is the next step).

wilmott <- brick("...air.mon.mean.v501.nc")
ls()
list2env(setNames(unstack(wilmott), names(wilmott)), .GlobalEnv)
ls()
lights <- raster("...F182013.v4/F182013.v4c_web.stable_lights.avg_vis.tif")

Both rasters have the following properties..

wilmott
class      : RasterBrick 
dimensions : 360, 720, 259200, 1416  (nrow, ncol, ncell, nlayers)
resolution : 0.5, 0.5  (x, y)
extent     : 0, 360, -90, 90  (xmin, xmax, ymin, ymax)
crs        : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
source     : ...air.mon.mean.v501.nc 
names      : X1900.01.01, X1900.02.01, X1900.03.01, X1900.04.01, X1900.05.01, X1900.06.01, X1900.07.01, X1900.08.01, X1900.09.01, X1900.10.01, X1900.11.01, X1900.12.01, X1901.01.01, X1901.02.01, X1901.03.01, ... 
Date/time  : 1900-01-01, 2017-12-01 (min, max)
varname    : air 
lights
class      : RasterLayer 
dimensions : 16801, 43201, 725820001  (nrow, ncol, ncell)
resolution : 0.008333333, 0.008333333  (x, y)
extent     : -180.0042, 180.0042, -65.00417, 75.00417  (xmin, xmax, ymin, ymax)
crs        : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0 
source     : ...F182013.v4/F182013.v4c_web.stable_lights.avg_vis.tif 
names      : F182013.v4c_web.stable_lights.avg_vis 
values     : 0, 255  (min, max)

and look like this Plot Temperature Data Plot Lights Data

As you can see both plots have a different extent. And now comes the problem. If I project (or resample (or playing around with the CRS, extents etc.)) my climate raster, I am able to match their CRS, extents etc. but a lot of data gets thrown away. See in the plot below:

wilmott.new.1901.01 = projectRaster(wilmott1901.01, lights) 
> wilmott.new.1901.01
class      : RasterLayer 
dimensions : 16801, 43201, 725820001  (nrow, ncol, ncell)
resolution : 0.008333333, 0.008333333  (x, y)
extent     : -180.0042, 180.0042, -65.00417, 75.00417  (xmin, xmax, ymin, ymax)
crs        : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0 
source     : /private/var/folders/dr/4gjhd01x5sd4wf965z0v11200000gp/T/RtmprJJ26D/raster/r_tmp_2021-01-12_164836_576_12591.grd 
names      : X1901.01.01 
values     : -56.09999, 35.79999  (min, max)

Plot of climate data projected to light data
How can I solve this? Any ideas? And how would I proceed to retrieve only the values for all grids in both sources in a new data frame? With rastertoPoints()or extract()?

Hope you get what I mean. Thanks a lot!


与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
Welcome To Ask or Share your Answers For Others

1 Answer

0 votes
by (71.8m points)
等待大神答复

与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
Welcome to OStack Knowledge Sharing Community for programmer and developer-Open, Learning and Share
Click Here to Ask a Question

...