rm(list=ls())
library(stargazer)
library(data.table)
library(ggplot2)
library(lfe)
library(rgdal)
library(rgeos)
library(gridExtra)
library(scales)
library(zoo)
library(foreign)
library(fst)
library(stringr)
library(ggrepel)
library(lubridate)
library(tidycensus)
library(arrow)
library(dplyr)
library(purrr)
library(jsonlite)
library(quantmod)
library(tidyverse)
library(R.matlab)
# library(directlabels)
fips <- fread("fips.csv")
names(fips) <- "fips"
states <- floor(fips/1000)
statefips <- fread("C:/Users/dratnadiwakara2/Documents/OneDrive - Louisiana State University/Raw Data/Crosswalk Files/us-state-ansi-fips.csv")
statefips <- statefips[,c("statefips","statecode")]
states <- merge(states,statefips,by.x="fips",by.y="statefips")
statenames <- c('AL','AK','AZ','AR','CA','CO','CT','DE','DC','FL','GA','ID','IL','IN','IA','KS','KY','LA','ME','MD','MA','MI','MN','MS','MO','MT','NE','NV','NH','NJ','NM','NY','NC','ND','OH','OK','OR','PA','RI','SC','SD','TN','TX','UT','VT','VA','WA','WV','WI','WY')
bmat <- fread("bmaterra.csv")
names(bmat) <- paste0("week",c("01","02","03","04","05","06","07","08","09","10","11","12"))
statefe <- bmat[1:50,]
statefe <- cbind(statenames,statefe)
names(statefe) <- c("State", paste0("week",c("01","02","03","04","05","06","07","08","09","10","11","12")))
statefe <- melt(statefe,id.vars = "State")
exovars <- bmat[51:55,]
exovars[,var:=c("pop","child","hs","white","inc")]
exovars <- melt(exovars,id.var="var")
ggplot(exovars,aes(x=variable,y=value,color=var,group=var))+geom_line()+theme_minimal()+theme(legend.position="bottom",legend.title = element_blank())+ xlab("")+ylab("")+labs(title="Exogenous Variables (Dep. Var: Pct Home)",subtitle="SDEM Estimates")

# ggplot(statefe,aes(x=variable,y=value)) + geom_text(aes(label=State),hju st=0, vjust=0,size=2)+theme_minimal()
ggplot(statefe,aes(x=variable,y=value))+geom_boxplot()+ geom_text(aes(label=State),hjust=1.5, vjust=0,size=1.5,check_overlap = TRUE)+theme_minimal()+ xlab("")+ylab("")+labs(title="State Fixed Effects (Dep. Var: Pct Home)",subtitle="SDEM Estimates")

fips <- fread("fips.csv")
names(fips) <- "fips"
xsub <- fread("xsub.csv")
names(xsub) <- c("pop","child","hs","white","inc")
pcthome <- fread("pcthome.csv")
pcthome <- log(pcthome)
pcthome <- pcthome - pcthome$V1
names(pcthome) <- paste0("pcthome_week",c("00","01","02","03","04","05","06","07","08","09","10","11","12"))
df <- cbind(pcthome,fips,xsub)
df[,state:=floor(fips/1000)]
df <- merge(df,statefips,by.x="state",by.y="statefips")
statefe <- NULL
coefs <- NULL
for(w in c("01","02","03","04","05","06","07","08","09","10","11","12") ) {
reg <- felm(as.formula(paste0("pcthome_week",w,"~pop+child+hs+white+inc|statecode")),data=df)
# t <- summary(reg)
# rsq <-c(rsq,t$adj.r.squared)
# print(t$adj.r.squared)
temp <- data.frame(getfe(reg)[,c("effect","idx")])
temp$effect <- as.numeric(temp$effect)
temp$idx <- as.character(temp$idx)
names(temp) <- c(paste0("pcthome_week",w),"state")
temp2 <- data.frame(reg$coefficients)
temp2['var'] <- row.names(temp2)
row.names(temp2) <- NULL
if(is.null(statefe)) {
statefe = temp
coefs <- temp2
} else{
statefe <- merge(statefe,temp,by="state")
coefs <- merge(coefs,temp2,by="var")
}
}
statefe <- melt(statefe,id.vars = "state")
coefs <- melt(coefs,id.vars = "var")
# ggplot(statefe,aes(x=variable,y=value)) + geom_text(aes(label=state),hjust=0, vjust=0,size=2)+theme_minimal()
ggplot(coefs,aes(x=substr(variable,13,14),y=value,color=var,group=var))+geom_line()+theme_minimal()+theme(legend.position="bottom",legend.title = element_blank())+ xlab("Week")+ylab("")+labs(title="Exogenous Variables (Dep. Var: Pct Home)",subtitle="OLS Estimates")

ggplot(statefe,aes(x=substr(variable,13,14),y=value))+geom_boxplot()+ geom_text(aes(label=state),hjust=1.5, vjust=0,size=1.5,check_overlap = TRUE)+theme_minimal()+theme(legend.position="bottom",legend.title = element_blank())+ xlab("Week")+ylab("")+labs(title="State Fixed Effects (Dep. Var: Pct Home)",subtitle="OLS Estimates")

# .printtable(reg)
zip_county <- fread("C:/Users/dratnadiwakara2/Documents/OneDrive - Louisiana State University/Raw Data/Crosswalk Files/ZIP_COUNTY_092016.csv")
zip_county<- zip_county[order(zip_county$ZIP,-zip_county$RES_RATIO)]
zip_county <- zip_county[!duplicated(zip_county$ZIP),c("ZIP","COUNTY")]
place_categories <- list()
place_categories[['grocery']] = c("General Merchandise Stores, including Warehouse Clubs and Supercenters","Grocery Stores")
gasolineother = c("Gasoline Stations","Health and Personal Care Stores")
place_categories[['gasolineother']] = c("Gasoline Stations","Health and Personal Care Stores")
place_categories[['restaurants']] =c("Restaurants and Other Eating Places")
place_categories[['entertainment1']] = c("Museums, Historical Sites, and Similar Institutions","Gambling Industries","Motion Picture and Video Industries","Spectator Sports")#,"Department Stores"
place_categories[['entertainment2']] = c("Sporting Goods, Hobby, and Musical Instrument Stores","Amusement Parks and Arcades")
place_categories[['healthservices']] = c("General Medical and Surgical Hospitals","Offices of Other Health Practitioners")
place_categories[['airports']] = c("Support Activities for Air Transportation")
core = read_fst(
"C:/Users/dratnadiwakara2/Documents/OneDrive - Louisiana State University/Raw Data/COVID19/sg/core.fst",
columns = c("safegraph_place_id","location_name","top_category","naics_code","latitude","longitude","region","postal_code","iso_country_code"),
as.data.table = TRUE)
files <- list.files(path = "C:/Users/dratnadiwakara2/Documents/OneDrive - Louisiana State University/Raw Data/COVID19/sg/weekly",pattern = '\\.fst',full.names = TRUE)
sg = lapply(files, read_fst, columns = c("safegraph_place_id","location_name","postal_code","naics_code","date_range_start","raw_visit_counts","raw_visitor_counts",paste0("visitsonday",1:7),"android","ios"),
as.data.table = TRUE)
sg <- do.call(rbind , sg)
sg = merge(sg,core[,c("safegraph_place_id","top_category")],by="safegraph_place_id")
sg = merge(sg,zip_county,by.x="postal_code",by.y="ZIP",all.x=TRUE)
sg[,week:=as.numeric(format(date_range_start,format="%Y%m%d"))]
sg[,place_category:=NA]
sg[,category:=ifelse(top_category %in% place_categories[['grocery']],"grocery",
ifelse(top_category %in% place_categories[['restaurants']],"restaurants",
ifelse(top_category %in% place_categories[['entertainment1']],"entertainment",NA )))]
wkstart <- c("3_01","3_08","3_15","3_22","3_29","4_05","4_12","4_19","4_26")
category_county_week = sg[,.(v1=sum(visitsonday1,na.rm=TRUE),
v2=sum(visitsonday2,na.rm=TRUE),
v3=sum(visitsonday3,na.rm=TRUE),
v4=sum(visitsonday4,na.rm=TRUE),
v5=sum(visitsonday5,na.rm=TRUE),
v6=sum(visitsonday6,na.rm=TRUE),
v7=sum(visitsonday7,na.rm=TRUE)),
by=.(category,week,COUNTY)]
category_county_week <- category_county_week[!is.na(category) & !is.na(COUNTY)]
category_county_week[,visits:=v1+v2+v3+v4+v5+v6+v7]
category_county_week[,c("v1","v2","v3","v4","v5","v6","v7"):=list(NULL)]
grocery <- dcast(category_county_week[category=="grocery",c("COUNTY","week","visits")],COUNTY~week,value.var = "visits")
grocery[,2:ncol(grocery)] <- log(1+grocery[,2:ncol(grocery)] )
grocery[,2:ncol(grocery)] <- grocery[,2:ncol(grocery)] - grocery$`20200301`
names(grocery) <- c("fips",paste0("grocery_week",wkstart))
df <- merge(df,grocery,by="fips",all.x=T)
restaurants <- dcast(category_county_week[category=="restaurants",c("COUNTY","week","visits")],COUNTY~week,value.var = "visits")
restaurants[,2:ncol(restaurants)] <- log(1+restaurants[,2:ncol(restaurants)] )
restaurants[,2:ncol(restaurants)] <- restaurants[,2:ncol(restaurants)] - restaurants$`20200301`
names(restaurants) <- c("fips",paste0("restaur_week",wkstart))
df <- merge(df,restaurants,by="fips",all.x=T)
entertainment <- dcast(category_county_week[category=="entertainment",c("COUNTY","week","visits")],COUNTY~week,value.var = "visits")
entertainment[,2:ncol(entertainment)] <- log(1+entertainment[,2:ncol(entertainment)] )
entertainment[,2:ncol(entertainment)] <- entertainment[,2:ncol(entertainment)] - entertainment$`20200301`
names(entertainment) <- c("fips",paste0("enterta_week",wkstart))
df <- merge(df,entertainment,by="fips",all.x=T)
cats <- c("grocery","restaur","enterta")
titles <- c("Grocery Visits","Restaurant Visits","Entertainment visits")
i=1
for(cat in cats) {
statefe <- NULL
coefs <- NULL
rsq <- NULL
gr1<- list()
gr2<-list()
print(titles[i])
for(w in wkstart[2:length(wkstart)] ) {
reg <- felm(as.formula(paste0(cat,"_week",w,"~pop+child+hs+white+inc|statecode")),data=df)
t <- summary(reg)
rsq <-c(rsq,t$adj.r.squared)
temp <- data.frame(getfe(reg)[,c("effect","idx")])
temp$effect <- as.numeric(temp$effect)
temp$idx <- as.character(temp$idx)
names(temp) <- c(paste0("grocery_week",w),"state")
temp2 <- data.frame(reg$coefficients)
temp2['var'] <- row.names(temp2)
row.names(temp2) <- NULL
if(is.null(statefe)) {
statefe = temp
coefs <- temp2
} else{
statefe <- merge(statefe,temp,by="state")
coefs <- merge(coefs,temp2,by="var")
}
}
statefe <- melt(statefe,id.vars = "state")
coefs <- melt(coefs,id.vars = "var")
# ggplot(statefe,aes(x=variable,y=value)) + geom_text(aes(label=state),hjust=0, vjust=0,size=2)+theme_minimal()
print(gr1[[i]] <- ggplot(coefs,aes(x=substr(variable,13,16),y=value,color=var,group=var))+geom_line(size=1)+geom_point(size=2)+theme_minimal()+theme(legend.position="bottom",legend.title = element_blank())+ xlab("Week")+ylab("")+labs(title=paste0("Exogenous Variables (Dep. Var: ",titles[i],")"),subtitle="OLS Estimates"))
print(gr2[[i]] <- ggplot(statefe,aes(x=substr(variable,13,16),y=value))+geom_boxplot()+ geom_text(aes(label=state),hjust=1.5, vjust=0,size=2,check_overlap = TRUE)+theme_minimal()+theme(legend.position="bottom",legend.title = element_blank())+ xlab("Week")+ylab("")+labs(title=paste0("State Fixed Effects (Dep. Var: ",titles[i],")"),subtitle="OLS Estimates"))
i=i+1
}
## [1] "Grocery Visits"


## [1] "Restaurant Visits"


## [1] "Entertainment visits"

