library(tidycensus)

census_api_key("06b232797e7854aad802d7c7d8673c337fe1b29a")


acs <- get_acs(geography = "county",variables = c("B06012_002","B01001_001"),year = 2018)

acs<-data.table(acs)
acs[,c("NAME","moe"):=list(NULL)]

acs <- dcast(acs,GEOID~variable)
acs[,pctbelowpoverty:=B06012_002/B01001_001]
acs[,GEOID:=as.numeric(GEOID)]
zipcity <- fread("C:/Users/dratnadiwakara2/Documents/OneDrive - Louisiana State University/Raw Data/Crosswalk Files/zip_code_mapping.csv")

tbill3month <- fread("C:/Users/dratnadiwakara2/Documents/OneDrive - Louisiana State University/Projects/COVID19/Muni/WTB3MS.csv")
tbill3month[,DATE:=as.Date(DATE)]
tbill3month[,week:=cut(DATE, "week")] 
tbill3month[,week:=as.character(week)]

zhvi <- readRDS("C:/Users/dratnadiwakara2/Documents/OneDrive - Louisiana State University/Raw Data/Zillow Research Data/Zip_Zhvi_SingleFamilyResidence.rds")
zhvi <- data.table(zhvi)
zhvi <- zhvi[month>="2019-01-01"]
zhvi <- zhvi[,.(zhvi=mean(zhvi)),by=(zipcode)]
zhvi[,zipcode:=as.integer(zipcode)]

zri <- readRDS("C:/Users/dratnadiwakara2/Documents/OneDrive - Louisiana State University/Raw Data/Zillow Research Data/Zip_Zri_SingleFamilyResidenceRental.rds")
zri <- data.table(zri)
zri <- zri[month>="2019-01-01"]
zri <- zri[,.(zri=mean(mean_rent)),by=(zipcode)]
zri[,zipcode:=as.integer(zipcode)]

zip_county <- fread("C:/Users/dratnadiwakara2/Documents/OneDrive - Louisiana State University/Raw Data/Crosswalk Files/ZIP_COUNTY_092016.csv")
setorder(zip_county,ZIP,-RES_RATIO)
zip_county <- zip_county[!duplicated(zip_county$ZIP),c("ZIP","COUNTY")]
names(zip_county) <- c("zip","fips")


citydata <- fread("C:/Users/dratnadiwakara2/Documents/OneDrive - Louisiana State University/Raw Data/COVID19/cities_by_income_2.csv")
citydata[,highincome:=ifelse(income>=median(income),1,0)]

citydata <- merge(citydata,zhvi,by.x="lzip",by.y="zipcode",all.x=T)
citydata <- merge(citydata,zri,by.x="lzip",by.y="zipcode",all.x=T)
citydata <- merge(citydata,zip_county,by.x="lzip",by.y="zip",all.x=T)

spatial <- R.matlab::readMat('C:/Users/dratnadiwakara2/Documents/OneDrive - Louisiana State University/Raw Data/COVID19/mat/spatial_data.mat')
fips <- unlist(spatial$FIPS, use.names=FALSE)
pop <- unlist(spatial$Pop, use.names=FALSE)
land <- unlist(spatial$Land,use.names = FALSE)
acs2 <- data.frame(fips=as.numeric(fips),countypop=pop,countyland=land)

citydata <- merge(citydata,acs2,by="fips",all.x=T)
citydata[,popdensity:=citydata$countypop/citydata$countyland]



revsource <- fread("C:/Users/dratnadiwakara2/Documents/OneDrive - Louisiana State University/Projects/COVID19/Muni/rev_source.csv")
revsource[,c("issuer","fyear"):=list(NULL)]
revsource <- revsource[Ticker!=""]

citydata <- merge(citydata,revsource,by="Ticker",all.x=T)

citydata[,proptax_rev:=PropTax_Rev/Total_Revenues]

load(file="C:/Users/dratnadiwakara2/Documents/OneDrive - Louisiana State University/Raw Data/Presidential Election/countypres_2000-2016.RData")
election <- x
rm(x)
election <- data.table(election)
election <- election[election$year %in% c(2016)& election$party=="republican"]
election <- election[,c("year","FIPS","candidatevotes","totalvotes")]
election[,rep_pct:=election$candidatevotes/election$totalvotes]
election <- election[,c("FIPS","rep_pct")]

citydata <- merge(citydata,election,by.x="fips",by.y="FIPS",all.x = T)

citydata <- merge(citydata,acs,by.x="fips",by.y="GEOID")
muni <- fread("C:/Users/dratnadiwakara2/Documents/OneDrive - Louisiana State University/Projects/COVID19/Muni/large_small_cities_by_state.csv")
muni[,Date:=as.Date(Date,format="%m/%d/%Y")]

muni[,week:=cut(Date, "week")] 
muni[,weekno:= as.factor(week)]
muni[,weekno:= as.numeric(weekno)-39]
muni[,week:=as.character(week)] 

muni <- merge(muni,tbill3month,by="week")

muni[,excessyld:=ytm_avg-WTB3MS]


citydata[,c("V1","statefips","population"):=list(NULL)]
muni <- merge(muni,citydata,by=c("Ticker"))


munifinancials <- fread("C:/Users/dratnadiwakara2/Documents/OneDrive - Louisiana State University/Projects/COVID19/Muni/city_financials.csv")

muni <- merge(muni,munifinancials,by.x="Issuer",by.y="Issuer",all.x=T)
muni[,cash_assets:=cash/totassets]
muni[,cash_opex:=cash/totopexp]
muni[,proptax_lib:=pptytax/totliabilities]
muni[,proptax_opex:=pptytax/totopexp]
muni[,facweekno:=as.factor(weekno)]

1 YTM - 3 Month T-Bill by Maturity

1.1 Mean

munisum <- muni[maturity_id<6,.(yld=mean(excessyld,na.rm=T),yldsd=sd(excessyld,na.rm=T),yldq1=quantile(excessyld,0.25,na.rm=T),yldq3=quantile(excessyld,0.75,na.rm=T)),by=.(weekno,maturity_id)]

ggplot(munisum,aes(x=weekno,y=yld,color=factor(maturity_id),group=factor(maturity_id)))+geom_line()+theme_minimal()+theme(legend.position="bottom",axis.text.x=element_text(angle=90))

1.2 Standard deviation

ggplot(munisum,aes(x=weekno,y=yldsd,color=factor(maturity_id),group=factor(maturity_id)))+geom_line()+theme_minimal()+theme(legend.position="bottom",axis.text.x=element_text(angle=90))

1.3 Interquartile range

ggplot(munisum,aes(x=weekno,y=(yldq3-yldq1),color=factor(maturity_id),group=factor(maturity_id)))+geom_line()+theme_minimal()+theme(legend.position="bottom",axis.text.x=element_text(angle=90))

2 Regression

\[ excessyld_{i,w} = \Sigma_w \beta_{1w} \times log(\frac{proptax}{revenue}) \times w + \Sigma_w \beta_{2w} \times log(republican pct 2016) \times w \\+ \Sigma_w \beta_{3w} \times log(popdensity) \times w + \Sigma_w \beta_{4w} \times log(\frac{cash}{assets}) \times w + \Sigma_w \beta_{5w} \times log(Houseprice) \times w \\+ \Sigma_w \beta_{6w} \times log(\frac{Populationbelowpoverty}{Population}) \times w + \Sigma_w \beta_{7w} \times log(medianincome) \times w + City FE + Week FE\]

r <- list()
r[[1]] <- felm(excessyld~log(proptax_rev)*facweekno+log(rep_pct)*facweekno+log(popdensity)*facweekno+factor(weekno)*log(cash_assets)+log(zhvi)*factor(weekno)+log(income.x)*factor(weekno)|Ticker+weekno|0|Ticker,data=muni[weekno>= (-20) & maturity_id==2])
r[[2]] <- felm(excessyld~log(cash_assets)*facweekno+log(rep_pct)*facweekno+log(popdensity)*facweekno+log(proptax_rev)*factor(weekno)+log(zhvi)*factor(weekno)+log(income.x)*factor(weekno)|Ticker+weekno|0|Ticker,data=muni[weekno>= (-20) & maturity_id==2])
r[[3]] <- felm(excessyld~log(zhvi)*facweekno+log(rep_pct)*facweekno+log(popdensity)*facweekno+log(proptax_rev)*factor(weekno)+log(cash_assets)*factor(weekno)+log(income.x)*factor(weekno)|Ticker+weekno|0|Ticker,data=muni[weekno>= (-20) & maturity_id==2])
r[[4]] <- felm(excessyld~log(income.x)*facweekno+log(rep_pct)*facweekno+log(popdensity)*facweekno+log(proptax_rev)*factor(weekno)+log(cash_assets)*factor(weekno)+log(zhvi)*factor(weekno)|Ticker+weekno|0|Ticker,data=muni[weekno>= (-20) & maturity_id==2])
r[[5]] <- felm(excessyld~log(popdensity)*facweekno+log(rep_pct)*facweekno++log(income.x)*facweekno+log(proptax_rev)*factor(weekno)+log(cash_assets)*factor(weekno)+log(zhvi)*factor(weekno)|Ticker+weekno|0|Ticker,data=muni[weekno>= (-20) & maturity_id==2])
r[[6]] <- felm(excessyld~log(rep_pct)*facweekno+log(popdensity)*facweekno+log(income.x)*facweekno+log(proptax_rev)*factor(weekno)+log(cash_assets)*factor(weekno)+log(zhvi)*factor(weekno)|Ticker+weekno|0|Ticker,data=muni[weekno>= (-20) & maturity_id==2])
r[[7]] <- felm(excessyld~log(pctbelowpoverty)*facweekno+log(rep_pct)*facweekno+log(popdensity)*facweekno+log(income.x)*facweekno+log(proptax_rev)*factor(weekno)+log(cash_assets)*factor(weekno)+log(zhvi)*factor(weekno)|Ticker+weekno|0|Ticker,data=muni[weekno>= (-20) & maturity_id==2])

# .printtable(r)
r5 <- list()
r5[[1]] <- felm(excessyld~log(proptax_rev)*facweekno+log(rep_pct)*facweekno+log(popdensity)*facweekno+factor(weekno)*log(cash_assets)+log(zhvi)*factor(weekno)+log(income.x)*factor(weekno)|Ticker+weekno|0|Ticker,data=muni[weekno>= (-20) & maturity_id==5])
r5[[2]] <- felm(excessyld~log(cash_assets)*facweekno+log(rep_pct)*facweekno+log(popdensity)*facweekno+log(proptax_rev)*factor(weekno)+log(zhvi)*factor(weekno)+log(income.x)*factor(weekno)|Ticker+weekno|0|Ticker,data=muni[weekno>= (-20) & maturity_id==5])
r5[[3]] <- felm(excessyld~log(zhvi)*facweekno+log(rep_pct)*facweekno+log(popdensity)*facweekno+log(proptax_rev)*factor(weekno)+log(cash_assets)*factor(weekno)+log(income.x)*factor(weekno)|Ticker+weekno|0|Ticker,data=muni[weekno>= (-20) & maturity_id==5])
r5[[4]] <- felm(excessyld~log(income.x)*facweekno+log(rep_pct)*facweekno+log(popdensity)*facweekno+log(proptax_rev)*factor(weekno)+log(cash_assets)*factor(weekno)+log(zhvi)*factor(weekno)|Ticker+weekno|0|Ticker,data=muni[weekno>= (-20) & maturity_id==5])
r5[[5]] <- felm(excessyld~log(popdensity)*facweekno+log(rep_pct)*facweekno++log(income.x)*facweekno+log(proptax_rev)*factor(weekno)+log(cash_assets)*factor(weekno)+log(zhvi)*factor(weekno)|Ticker+weekno|0|Ticker,data=muni[weekno>= (-20) & maturity_id==5])
r5[[6]] <- felm(excessyld~log(rep_pct)*facweekno+log(popdensity)*facweekno+log(income.x)*facweekno+log(proptax_rev)*factor(weekno)+log(cash_assets)*factor(weekno)+log(zhvi)*factor(weekno)|Ticker+weekno|0|Ticker,data=muni[weekno>= (-20) & maturity_id==5])
r5[[7]] <- felm(excessyld~log(pctbelowpoverty)*facweekno+log(rep_pct)*facweekno+log(popdensity)*facweekno+log(income.x)*facweekno+log(proptax_rev)*factor(weekno)+log(cash_assets)*factor(weekno)+log(zhvi)*factor(weekno)|Ticker+weekno|0|Ticker,data=muni[weekno>= (-20) & maturity_id==5])
# .printtable(r)

2.1 Property Tax/Revenue

.coef_plot_2reg_line(r[[1]],"Maturity ID=2",r5[[1]],"Maturity ID=5","log(proptax_rev):facweekno",-20)

2.2 Cash/Assets

.coef_plot_2reg_line(r[[2]],"Maturity ID=2",r5[[2]],"Maturity ID=5","log(cash_assets):facweekno",-20)

2.3 House Price

.coef_plot_2reg_line(r[[3]],"Maturity ID=2",r5[[3]],"Maturity ID=5","log(zhvi):facweekno",-20)

2.4 Income

.coef_plot_2reg_line(r[[4]],"Maturity ID=2",r5[[4]],"Maturity ID=5","log(income.x):facweekno",-20)

2.5 Population Density

.coef_plot_2reg_line(r[[5]],"Maturity ID=2",r5[[5]],"Maturity ID=5","log(popdensity):facweekno",-20)

2.6 Pct of Republican 2016

.coef_plot_2reg_line(r[[6]],"Maturity ID=2",r5[[6]],"Maturity ID=5","log(rep_pct):facweekno",-20)

2.7 Pct Below Poverty

.coef_plot_2reg_line(r[[7]],"Maturity ID=2",r5[[7]],"Maturity ID=5","log(pctbelowpoverty):facweekno",-20)