setwd("~/Downloads")
library(ipumsr)
ddi <- read_ipums_ddi("usa_00008.xml")
data <- read_ipums_micro(ddi)
## Use of data from IPUMS USA is subject to conditions including that users should
## cite the data appropriately. Use command `ipums_conditions()` for more details.
data<-haven::zap_labels(data) #necessary to avoid problems with "labelled" data class
names(data)<-tolower(names(data))
options(tigris_class = "sf")
pumas<-pumas(state = "CA",
             year = 2019,
             cb = T)
data$pwt <- data$perwt
data$hwt <- data$hhwt
data$migration <- Recode(data$migrate1, recodes = "1=1; 0=NA; else=0")
des<-svydesign(ids = ~cluster,
               strata = ~ strata,
               weights = ~pwt,
               data = data)
#test<-svytable(~I(cit=="US born")+puma+sex, design=des )

puma_est_migration<-svyby(formula = ~migration,
                    by = ~puma,
                    design = des,
                    FUN=svymean,
                    na.rm = TRUE )
head(puma_est_migration)
##     puma migration          se
## 101  101 0.7442398 0.009366589
## 102  102 0.8259565 0.006778309
## 103  103 0.8749223 0.006679989
## 104  104 0.8904980 0.008340314
## 105  105 0.8888690 0.006785526
## 106  106 0.9227006 0.006058942

join to geography

pumas$puma<-as.numeric(pumas$PUMACE10)

geo1<-left_join(pumas, puma_est_migration, by=c("puma"= "puma"))
#head(geo1)
tmap_mode("plot")
## tmap mode set to plotting
tm_basemap("OpenStreetMap.Mapnik")+
  tm_shape(geo1)+
  tm_polygons("migration",
              style="kmeans",
              title=c("%  in the same house"),
              palette="Blues",
              n=8,
              legend.hist = TRUE) +
  tm_layout(legend.outside = TRUE,
            title = "% in the same house last year  \n California PUMAs 2015-2019",
            title.size =1.5,
            legend.frame = TRUE,
            ) +  tm_compass(position = c("right","top")) + tm_format("World",
            legend.position =  c("left", "bottom"),
            main.title.position =c("center")) + tm_scale_bar(position = c("left","bottom"))

tmap_mode("view")
## tmap mode set to interactive viewing
tm_basemap("OpenStreetMap.Mapnik")+
  tm_shape(geo1)+
  tm_polygons("migration",
              style="kmeans",
              title=c("%  in the same house"),
              palette="Blues",
              n=8,
              legend.hist = TRUE) +
  tm_layout(legend.outside = TRUE,
            title = "% in the same house last year  \n California PUMAs 2015-2019",
            title.size =1.5,
            legend.frame = TRUE,
            ) +  tm_compass(position = c("right","top")) + tm_format("World",
            legend.position =  c("left", "bottom"),
            main.title.position =c("center")) + tm_scale_bar(position = c("left","bottom"))
## Compass not supported in view mode.
## legend.postion is used for plot mode. Use view.legend.position in tm_view to set the legend position in view mode.