#I am trying to use GJRM with a Markov Random Field as explained here:
#https://www.rdocumentation.org/packages/GJRM/versions/0.2-6.5/topics/hiv
#I was able to replicate the study. But when I try using another shape file to create polygons, "GJRM"
#gives the following error("Error in if (all.equal(sort(a.name), sort(levels(k))) != TRUE) stop("mismatch #between #nb/polys supplied area names and data area names") : the condition has length > 1).
#However, when I test the "all. equal function" to the polygon names and the levels of the ID variable, they are #balanced.
#I don't know if a bug causes it or if I am making mistakes.
#The maintainer is not responding
#Example using fake data
library(tidyverse)
library(sf) #Create spatial object with simple features
library(GJRM) #Generalised Joint Regression Modelling
library(spdep) #Create Spatial Polygon Dataframe from sf
library(reshape2) #Data creation
#Create fake data with Zip Code Tabulated Areas (ZCTAs) from the year 2020
ZCTA5CE20 <- c("70001", "70002", "70003", "70005", "70006", "70030", "70031", "70032",
"70036", "70037", "70038", "70039", "70040", "70041")
year <- c("2019","2020", "2021")
df <- expand.grid(ZCTA5CE20=ZCTA5CE20, year=unique(year))
n <- nrow(df)
df$x <- runif(n = 14, min = 0, max = 1)
df$y <- df$x*runif(n = 14, min = 0, max = 1)
df <- df |> mutate(y1=ifelse(y>=0.5, 1, 0),
y2=ifelse(y==1 & y>=0.8, 1,0))
glimpse(df)
#########Create an SF Object with trigris package and filter
#the ZCTAs in the data
list_ZCTAS <- unique(df$ZCTA)
options(tigris_use_cache = TRUE)
spf_ZCTAS <- tigris::zctas(year=2020, class="sf") %>%
dplyr::filter(ZCTA5CE20 %in% list_ZCTAS ) %>% arrange(ZCTA5CE20)
#Add variables recognized by the "setup.polys" in the JGRM package
n=nrow(spf_ZCTAS
spf_ZCTAS$ID_0=244
spf_ZCTAS$ISO = "USA"
spf_ZCTAS$NAME_0 = "United States of America"
spf_ZCTAS$ID_1 =seq(1,n,1)
spf_ZCTAS$NAME_1=spf_ZCTAS$ZCTA5CE20
#Create a SpatialPolygonsDataFrame
spf_ZCTAS_sp <- as(spf_ZCTAS, "Spatial")
class(spf_ZCTAS_sp)
#Create polygons with vertices
row.names(spf_ZCTAS_sp)<- spf_ZCTAS_sp@data$NAME_1
spf_ZCTAS_polys <- GJRM::polys.setup(spf_ZCTAS_sp)
#ZCTAs as factor and check if equal
df$ZCTA5CE20 <- as.factor(df$ZCTA5CE20)
all.equal(sort(names(spf_ZCTAS_polys$polys)),
sort(levels(df$ZCTA5CE20)))
#Polygon list
xt= list(polys=spf_ZCTAS_polys$polys)
#GJRM Model
eq1 <- y1 ~ year + y + x
eq2 <- y2 ~ year + x
theta.eq <- ~ s(ZCTA5CE20, bs = "mrf", k=2, xt =xt )
f <- list(eq1, eq2, theta.eq)
bss <- gjrm(f, data = df, copula = "J90", model = "BSS",
margins = c("probit", "probit"))
#Error in if (all.equal(sort(a.name), sort(levels(k))) != TRUE) stop("mismatch between nb/polys supplied area #names and data area names") : the condition has length > 1
#But the test above was "TRUE", any help will be appreciated
```r