[R-sig-Geo] add spatraster layers by matching a dataframe column and make spatio-temporal rast

Julien Rodriguez ju||enrodr|guez@@rd@@pm @end|ng |rom gm@||@com
Mon May 13 18:35:11 CEST 2024


 Dear Marta,


I hope I've understood properly what you're expecting.

Provided your vector values are ordered properly, you can add it as a new
layer in an existing SpatRaster using "$" as for a data.frame.
The terra::merge method may be used for different SpatRasters or between
SpatVector and data.frames but you could prepare new layers from the
original raster values using base::merge with your "ID" column.
Regarding the temporal aspect, it would mean, as I understood it, that you
should have two data.frames containing the features to be included instead
of one (for years 2001 and 2002).
Surely not the most elegant way to do it below with base R functions but it
should work.
If anyone has a more appropriate solution, I'd be greatly interested too!


# runnable example
# Make an example raster:
require(terra)
set.seed(240513)

r <- terra::rast(ncols=10, nrows=10)
values(r) <- 1:100
names(r)="ID"
# add some missing values (to make it real)
r[c(2,43,60,94)]=NA
plot(r)

# Make an example data frame (note both have an ID col to match)
d <- data.frame(
  ID=100:1,
  k=c("aa","bb"),
  e=rnorm(100,2),
  year=sample(2001:2002, 100, replace = TRUE),
  date= as.Date(rep(c("2010-01-01", "2012-02-01",
                      "2010-03-01","2010-04-01"), l=100)))
head(d)

# the terra::time is used to describe different layers of the raster which
is not appropriated in this case (one layer = one feature for one year)
# Here the values to be included in the raster are associated to different
dates so the data frame has to be splitted by year
# A different year means actually a missing value for the pixel
# the raster time index to be then defined as:
ras.time.index <- 2001:2002

col2sel <- colnames(d)[ !colnames(d) %in% c("ID", "year")]

timed.d <- do.call(merge,
              lapply(ras.time.index , function(y){
                      new.d <- d
                      new.d [!d$year %in% y, col2sel] <- NA
                      colnames(new.d)[ colnames(new.d) %in% col2sel] <-
paste(col2sel, y, sep = "_")
                      new.d <- new.d[ order(new.d$ID), !colnames(new.d)
%in% "year"]
                      return(new.d)}))

head(timed.d)


# Extract the values from the raster and make a unique index to identify
the pixel values order from original raster
ras.vals <- terra::values(r, dataframe = TRUE)
ras.val.0 <- cbind(ras.vals, index = 1:nrow(ras.val))

# Merge with d (values with missing IDs will be lost from d) while keeping
every values from original raster
ras.val.1 <- merge(ras.val.0,  timed.d, by = "ID", all.x = TRUE)
summary(ras.val.1)

# Reorder features values properly
ras.val.1 <- ras.val.1[ order(ras.val.1$index), ]

# Add these new values to the initial raster iteratively
features <- colnames(timed.d)[ !colnames(timed.d) %in% colnames(ras.val.0)]
features

for(i in 1:length(features)){

  new.values <- ras.val.1 [, features[i]]
  if(inherits(new.values, "Date")){
    new.values <- as.character(new.values)
  }
  r$new.val <- new.values
  names(r)[ names(r) %in% "new.val"] <- features[i]

}

plot(r)

# At last, define the temporal component of different layers
# To do that, distinguishing the different features seems appropriated
# Example for feature "e" layers
ras.e <-  r[[substr(names(r), 1,2) %in% "e_"]]
terra::time(ras.e, tstep="years") <- ras.time.index

plot(ras.e)

I hope this might be useful

 Best regards

Julien Rodriguez







Le lun. 13 mai 2024 à 12:55, Marta Rufino <marta.m.rufino using gmail.com> a
écrit :

> Hi,
>
> I would like to add layers to a spatraster (terra), by matching (joining) a
> column in a data frame with the values of the raster layer (A).
> Then, I would like to make it a spatio-temporal terra spatraster and/or a
> stars (B).
>
> Any help will be strongly appreciated,
> Kind regards,
> M.
> ________________________________________
>  require(terra)
>
> # A)
>
> # runnable example
> # Make an example raster:
>   r = terra::rast(ncols=10, nrows=10)
>   values(r) <- 1:100
>   names(r)="ID"
>   # add some missing values (to make it real)
>   r[c(2,43,60,94)]=NA
>   plot(r)
>
>   # Make an example data frame (note both have an ID col to match)
>   d <- data.frame(
>     ID=100:1,
>     k=c("aa","bb"),
>     e=rnorm(100,2),
>     year=sample(2001:2002, 100, replace = TRUE),
>     date= as.Date(rep(c("2010-01-01", "2012-02-01",
> "2010-03-01","2010-04-01"), l=100)))
>   head(d)
>
> # My tries:
>
>   # 1. This works, but only for one col and it is not so nice:
>   r$year2 <- d$year[match(as.vector(r$ID), d$ID)]
>   plot(r)
>
>   # 2. try it using app function: I could not make it work
>   kk1 <- function(i,j=d) {
>       j[match(as.vector(i$ID), j$ID),]
>     }
>   app(r, fun = kk1)
>
>   # 3. try it using classify function
>   # works with numeric (not characters) but it is not perfect
>   terra::classify(r, d) # error
>   # two numeric cols: works, but add the second col as the first raster
> col. Not ok.
>   terra::classify(r, d[,c(1,3)])
>   # excluding the character col: gives an error
>   terra::classify(r, d[,-c(2)])
>
>   # 4. try it using merge function
>   terra::merge(r, d) %>% head() #not sure what he did and in which order :(
> but clearly not match/join working properly
>
>   # 5. only (not elegant) way I managed:
>   rast(
>     left_join(
>       as_tibble(r, xy=TRUE), d, by="ID"), type="xyz", crs = "EPSG:4326")
> %>% plot()
>
> # B) Make terra r raster a spatiotemporal cube
>
> # 1.
> #
>
> https://stackoverflow.com/questions/74079441/group-raster-files-by-week-month
> terra::time(r, tstep="years") <- d$year
> # gives an error...
> # Error: [time<-] length(value) != nlyr(x)
>
> # 2.
> #
>
> https://stackoverflow.com/questions/74934891/specifying-layers-when-using-rast-type-xyz-to-convert-data-frame-to-spatra
> # but I think this is not really a spatio-temporal rast...
> d2 <- left_join(as_tibble(r, xy=TRUE), d, by=c("ID"))
> w <- reshape(d2, timevar="year", idvar=c("x", "y"), direction="wide")
> x <- rast(w, type="xyz")
>
> # by year
> d3 <- split(d2[,-7], d2$year)
> r2 <- lapply(d3, \(i) rast(i, type="xyz"))
> (r2 <- rast(r2))
> time(r2) <- as.Date(names(r2))
> plot(r2)
>
> # plot with rasterVis
> require(rasterVis)
> levelplot(r2)
>
>
> # 2. using stars (aternative)
> # sources:
>
> https://stackoverflow.com/questions/77368957/how-to-transform-a-sf-table-to-a-stars-cube-with-starsst-as-stars-function
> # (
>
> https://gis.stackexchange.com/questions/353741/generate-a-stars-object-from-tabular-spatial-temporal-data
> ):
>
> s <- st_as_stars(as.data.frame(r, xy = TRUE, na.rm = TRUE), dims = c("x",
> "y", "year"), crs = 4326)
> s %>% plot(axes=TRUE)
> # but the from and to of the year are wrong
> # How to change to spatrast again?
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo using r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>

	[[alternative HTML version deleted]]



More information about the R-sig-Geo mailing list