Hi all,
I am trying to use PGOcc and spPGOcc to model marten occupancy over a section of National Forest land in Idaho. I have successfully fitted short simple runs of both models with a small number of occ.covs and det.covs to practice using the functions, and I can back transform to look at the occupancy estimates and detection estimates for my covariates. However, I am having trouble getting a prediction surface for occupancy using the covariates in my model.
I'm familiar with using terra::predict(raster, fitted.mod) to predict out across a raster surface. However, it seems that the model output from (sp)PGOcc require that you specify the coordinates to predict over as a df or matrix. I've gone back to try to get this to work on the non-spatial model first, shown here.
I can successfully generate an occupancy prediction by just giving the covariate values at the camera locations, which I'll attached the output figure for as well.
The problem seems to be the memory size generating the prediction at all coordinates of my intended raster. I'm interested if there is a workaround for this? Is there a way to predict the fitted model with a raster of interest as one can do with other fitted model objects?
Thanks!
Kristin
# package data. occ.covs contains several covariates extracted at camera locations and then scaled, including DEM.sc, which I am going to use in the model.
marten_data <- list( y=station_det_wk_bin,
occ.covs = occ.covs,
det.covs = list(effort=station_det_wk_eff,
Lure = station_metadata.df$Lure,
ScentFamily = station_metadata.df$ScentFamily,
Bait.Type = station_metadata.df$Bait.Type),
coords = station.coords3)
n.samples <- 3000
n.burn <- 1000
n.thin <- 2
n.chains <- 1
mod1 <- PGOcc(occ.formula = ~DEM.sc,
det.formula = ~Lure,
data = marten_data,
# inits = oven.inits,
n.samples = n.samples,
# priors = oven.priors,
n.omp.threads = 1,
verbose = TRUE,
n.report = 1000,
n.burn = n.burn,
n.thin = n.thin,
n.chains = n.chains)
summary(mod1)
summary(mod1)
Call:
PGOcc(occ.formula = ~DEM.sc, det.formula = ~Lure, data = marten_data, n.samples = n.samples, n.omp.threads = 1, verbose = TRUE, n.report = 1000,
n.burn = n.burn, n.thin = n.thin, n.chains = n.chains)
Samples per Chain: 3000
Burn-in: 1000
Thinning Rate: 2
Number of Chains: 1
Total Posterior Samples: 1000
Run Time (min): 0.0392
Occurrence (logit scale):
Mean SD 2.5% 50% 97.5% Rhat ESS
(Intercept) 0.7201 0.2149 0.3134 0.7187 1.1447 NA 661
DEM.sc 0.3549 0.2092 -0.0469 0.3542 0.7530 NA 861
Detection (logit scale):
Mean SD 2.5% 50% 97.5% Rhat ESS
(Intercept) -1.3776 0.1048 -1.5932 -1.3725 -1.1739 NA 1000
LureNone -0.0746 0.2981 -0.6727 -0.0535 0.4797 NA 604
LureSingleDose -0.0139 0.1409 -0.2761 -0.0233 0.2604 NA 887
# predict and plot data at camera locations
pred1.MOD2 <- predict(mod1, X.0=cbind(1,marten_data$occ.covs$DEM.sc))
plot_df.2 <- data.frame( x = marten_data$coords$long,
y = marten_data$coords$lat,
mean.psi = apply(pred1.MOD2$psi.0.samples,2,mean),
sd.psi = apply(pred1.MOD2$psi.0.samples, 2, sd))
dat.stars <- st_as_stars(plot_df.2, dims = c('x', 'y'))
head(dat.stars)
ggplot() +
geom_stars(data=dat.stars, aes(x, y, fill = mean.psi)) +
scale_fill_viridis_c(na.value = "transparent") +
labs(x = "Longitude", y = "Latitude", fill = "Predicted\nOccupancy")
* see attached fig*
# try to plot over area of interest; modify the original DEM layer into a df to use predict() and then scale it
coords.df <- crds(DEM, df=TRUE)
DEM.df <- as.data.frame(DEM)
elev.as.df <- data.frame(coords.df,DEM.df)
elev.pred <- (elev.as.df$Elev - mean(station.coords2$DEM$Elev) / sd(station.coords2$DEM$Elev))
X.0 <- cbind(1, elev.pred)
head(X.0)
head(X.0)
elev.pred
[1,] 1 1300.863
[2,] 1 1306.006
[3,] 1 1309.755
[4,] 1 1315.756
[5,] 1 1322.356
[6,] 1 1328.449
pred.MOD <- predict(mod1, X.0)
Error: cannot allocate vector of size 194.1 Gb