Rent Conley SE Regression Models
library(tidyverse)## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.2 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.0
## ✔ ggplot2 3.4.2 ✔ tibble 3.2.1
## ✔ lubridate 1.9.2 ✔ tidyr 1.3.0
## ✔ purrr 1.0.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorslibrary(maps)##
## Attaching package: 'maps'
##
## The following object is masked from 'package:purrr':
##
## maplibrary(mapdata)
library(lubridate)
library(readr)
library(dplyr)
library(stringr)
library(plm)##
## Attaching package: 'plm'
##
## The following objects are masked from 'package:dplyr':
##
## between, lag, leadlibrary(lmtest)## Loading required package: zoo
##
## Attaching package: 'zoo'
##
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numericlibrary(synthdid)
library(devtools)## Loading required package: usethislibrary(latex2exp)
library(MCPanel)
library(rngtools)
library(future)
library(doFuture)## Loading required package: foreach
##
## Attaching package: 'foreach'
##
## The following objects are masked from 'package:purrr':
##
## accumulate, whenlibrary(future.batchtools)## Loading required package: parallellylibrary(xtable)
library(tidyr)
library(tibble)
library(cowplot)##
## Attaching package: 'cowplot'
##
## The following object is masked from 'package:lubridate':
##
## stamplibrary(sf)## Linking to GEOS 3.11.0, GDAL 3.5.3, PROJ 9.1.0; sf_use_s2() is TRUElibrary(conleyreg)
library(gridExtra)##
## Attaching package: 'gridExtra'
##
## The following object is masked from 'package:dplyr':
##
## combineConverting Rent data to Balanced Data
The rent data was converted to a balanced data (i.e equal observations for all year), then it was converted to a panel matrix data to be used in SDID and SC Estimate.
RENT.merge<-read.csv("MSA-RENT merge All.csv")
RENT.merge$treat<-ifelse(RENT.merge$state=="LA",1,0)
colnames(RENT.merge)## [1] "FirstNAME" "AREANAME" "RENT_0"
## [4] "RENT_1" "RENT_2" "RENT_3"
## [7] "RENT_4" "YEAR" "NAME.x"
## [10] "state" "Unemployee.Rate" "Personal.Income"
## [13] "Resident.Population" "nyear" "MSA"
## [16] "MSA.Code" "Price" "Change"
## [19] "NAME.y" "treat"# Remove NA
RENT.merge<-RENT.merge %>%
group_by(MSA)%>%
filter(!any(is.na(Unemployee.Rate)))
#KEEP NEW ORLEAN
RENT.merge1<-RENT.merge
RENT.merge1$NS<-paste0(RENT.merge1$FirstNAME,RENT.merge1$state)
RENT.NEWLA<-RENT.merge1%>%
filter(NS=="newLA")
RENTUS<-RENT.merge1%>%
filter(state!="LA")
RENT.merge1<-RENT.merge1%>%
filter(NS!="newLA")
## WITH NEW ORLEAN ONLY
RENT.merge2<-rbind(RENT.NEWLA,RENTUS)
RENT.merge2$YEAR<-as.numeric(RENT.merge2$YEAR)
RENT.merge2$s.treat<-ifelse(RENT.merge2$state=="LA" &
RENT.merge2$YEAR>2005,1,0)
##BALLANCE BY NS
unq<-data.frame(NS=rep(unique(RENT.merge2$NS),15))
unq<-arrange(unq,NS)
unq$YEAR<-rep(2001:2015,length(unique(RENT.merge2$NS)))
Ballance.RENT<-left_join(unq,RENT.merge2,by=c("YEAR","NS"))
Ballance.RENT<-Ballance.RENT%>%
group_by(NS)%>%
filter(!any(is.na(RENT_0)))
Ballance.RENT$lg<-log(Ballance.RENT$RENT_0)
Ballance.RENT$NS<-factor(Ballance.RENT$NS)
Ballance.RENT<-data.frame(Ballance.RENT)
panelp<-panel.matrices(Ballance.RENT,unit = "NS",time = "YEAR",outcome = "lg",treatment = "s.treat")SDID Estimate
## SDID Estimate
r0<-synthdid_estimate(Y=panelp$Y,N0 = panelp$N0,T0 = panelp$T0)
plot(r0,overlay=1)+theme(legend.position = "bottom",
legend.direction = "horizontal",
legend.box.background = element_blank())+
ylab("RENT 0")+xlab("Year")
SC Estimate
## SC Estimate
sc0<-sc_estimate(Y=panelp$Y,N0 = panelp$N0,T0 = panelp$T0)
plot(sc0)+theme(legend.position = "bottom",
legend.direction = "horizontal",
legend.box.background = element_blank())+
ylab("RENT 0")+xlab("Year")
Unit and Time Weights Functions
The functions to be used in calculating the unit and time weights were created.
### Unit and Time Weights
mc_estimate = function(Y, N0, T0) {
N1=nrow(Y)-N0
T1=ncol(Y)-T0
W <- outer(c(rep(0,N0),rep(1,N1)),c(rep(0,T0),rep(1,T1)))
mc_pred <- mcnnm_cv(Y, 1-W, num_lam_L = 20)
mc_fit <- mc_pred$L + outer(mc_pred$u, mc_pred$v, '+')
mc_est <- sum(W*(Y-mc_fit))/sum(W)
mc_est
}
mc_placebo_se = function(Y, N0, T0, replications=200) {
N1 = nrow(Y) - N0
theta = function(ind) { mc_estimate(Y[ind,], length(ind)-N1, T0) }
sqrt((replications-1)/replications) * sd(replicate(replications, theta(sample(1:N0))))
}
difp_estimate = function(Y, N0, T0) {
synthdid_estimate(Y, N0, T0, weights=list(lambda=rep(1/T0, T0)), eta.omega=1e-6)
}
sc_estimate_reg = function(Y, N0, T0) {
sc_estimate(Y, N0, T0, eta.omega=((nrow(Y)-N0)*(ncol(Y)-T0))^(1/4))
}
difp_estimate_reg = function(Y, N0, T0) {
synthdid_estimate(Y, N0, T0, weights=list(lambda=rep(1/T0, T0)))
}
estimators <- list(did = did_estimate,
sc = sc_estimate,
sdid = synthdid_estimate,
difp = difp_estimate,
mc = mc_estimate,
sc_reg = sc_estimate_reg,
difp_reg = difp_estimate_reg)
##
estimates <- lapply(estimators, function(estimator) {
estimator(panelp$Y, panelp$N0, panelp$T0)
})Unit and Time Weights
unit.weights <- synthdid_controls(estimates[1:3], weight.type='omega', mass=1)
time.weights <- synthdid_controls(estimates[1:3], weight.type='lambda', mass=1)
unit.df <- as.data.frame(round(unit.weights[rev(1:nrow(unit.weights)), ], 3)) %>%
rownames_to_column(var = "NS")
time.df <- as.data.frame(round(time.weights, 3)) %>%
rownames_to_column(var = "NS")
head(unit.df, 10)## NS did sc sdid
## 1 abileneTX 0.006 0 0.000
## 2 akronOH 0.006 0 0.020
## 3 albanyGA 0.006 0 0.003
## 4 albuquerqueNM 0.006 0 0.000
## 5 altoonaPA 0.006 0 0.000
## 6 amarilloTX 0.006 0 0.008
## 7 anchorageAK 0.006 0 0.000
## 8 annistonAL 0.006 0 0.002
## 9 annMI 0.006 0 0.005
## 10 ashevilleNC 0.006 0 0.000Creating Data for Map plotting
MSA_code <- Ballance.RENT %>%
distinct(NS, MSA.Code)
# Adding MSA code to our unit weights data
unit.df1 <- inner_join(MSA_code, unit.df, by = "NS")
map_df <- unit.df1 %>%
mutate(sc_mp = ifelse(sc > 0, 1, 0),
sdid_mp = ifelse(sdid > 0, 1, 0))
# Loading the shape data for each MSA to be used in Mapping
land_shapefile <- "/Users/owad/Downloads/tl_2019_us_cbsa/tl_2019_us_cbsa.shp"
land_areas <- st_read(land_shapefile) %>%
rename(MSA.Code = CBSAFP) %>%
mutate(MSA.Code = as.numeric(MSA.Code))## Reading layer `tl_2019_us_cbsa' from data source
## `/Users/owad/Downloads/tl_2019_us_cbsa/tl_2019_us_cbsa.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 938 features and 12 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -178.4436 ymin: 17.83151 xmax: -65.42271 ymax: 65.45352
## Geodetic CRS: NAD83# Merging with weights data
msa_area <- inner_join(map_df, land_areas, by = "MSA.Code")
map_sdid <- msa_area %>%
filter(sdid_mp == 1)
map_sc <- msa_area %>%
filter(sc_mp == 1)Creating USA Map
usa_map <- map_data("state")
usa_plot <- ggplot() +
geom_polygon(data = usa_map, aes(x = long, y = lat, group = group), fill = "white", color = "black") +
labs(title = "Control MSAs") +
theme_map() +
theme(
plot.title = element_text(size = 16, hjust = 0.5),
axis.text = element_blank(),
axis.title = element_blank(),
panel.background = element_blank(),
panel.grid = element_blank(),
plot.margin = unit(c(1, 1, 1, 1), "cm")
)SDID Control MSAs
# Plot the MSA points on the USA map for SDID == 1
usa_plot +
geom_sf(data = map_sdid,
aes(geometry = geometry), fill = "lightblue", color = "black") +
theme(
plot.title = element_text(size = 18, hjust = 0.5),
plot.margin = unit(c(1, 1, 1, 1), "cm")
)
SC Control MSAs
usa_plot +
geom_sf(data = map_sc,
aes(geometry = geometry), fill = "lightgreen", color = "black") +
theme(
plot.title = element_text(size = 18, hjust = 0.5),
plot.margin = unit(c(1, 1, 1, 1), "cm")
)
Event Study for each Rent based on MSA with weight (sdid > 0)
## Event Study for MSAs with weight only and Louisiana MSAs
lat_long <- read.csv("rent data new.csv") %>%
distinct(X, Y, AREANAME)
RENT.merge <- RENT.merge %>%
merge(lat_long, by = "AREANAME")
# Lousiana MSAs
LA <- RENT.merge %>%
filter(state == "LA") %>%
distinct(MSA, MSA.Code)
# MSAs with weight (sdid > 0)
MSA_sdid <- map_sdid %>%
select(NAME, MSA.Code) %>%
rename(MSA = NAME)
## MSAs with weight (sdid > 0) and Louisiana MSAs
LA_sdid <- rbind(LA, MSA_sdid)
# Subsetting data to only MSAs with weight (sdid > 0) and Louisiana MSAs
# Data is then used for Event study
df <- inner_join(RENT.merge, LA_sdid, by = c("MSA", "MSA.Code"))
dloop <- data.frame(YEAR = unique(df$YEAR))
## YEAR DUMMY
for (i in 1:nrow(dloop)) {
df[paste0("d.", dloop$YEAR[i])] <- as.numeric(df$YEAR == unique(df$YEAR)[i])
}Rent 0
cse <- conleyreg(log(RENT_0) ~ Price + Unemployee.Rate +
Personal.Income + Resident.Population +
treat:(d.2001 + d.2002 + d.2003 +
d.2005 + d.2006 + d.2007 + d.2008 +
d.2009 + d.2010 + d.2011 + d.2012 +
d.2013 + d.2014 + d.2015) | MSA + YEAR,
data = df,
dist_cutoff = 300,
st_distance = TRUE,
lat = "Y",
lon = "X",
)## Checking spatial attributes## Estimating model## Warning: In fixest_env(fml = fml, data = data, weights = weig...:
## OpenMP not detected: cannot use 4 threads, single-threaded mode instead.## Estimating distance matrix and addressing spatial correlationcolnames(cse)[2] <- "Conley SE"
coef_df <- cse[,] %>%
as_tibble() %>%
mutate(variables = rownames(cse))
coef_df <- coef_df[5:18,]
event_df <- data.frame(time = 1:nrow(coef_df),
coef = coef_df$Estimate,
cse = coef_df$`Conley SE`)
event_df <- rbind(event_df, c("time" = 3.5, "coef" = 0, "cse" = 0))
event_df <- arrange(event_df, time)
# Calculate confidence intervals
event_df$ci_upper <- event_df$coef + 1.96 * event_df$cse
event_df$ci_lower <- event_df$coef - 1.96 * event_df$cse
# Extract dates from row names of coef_df
##dates <- as.Date(paste0(row.names(coef_df),"0101"), format = "treat:YEAR%Y%m%d")
dates<-seq.Date(from = as.Date("2001-01-01"),as.Date("2015-01-01"),by="year")
# Define start and end dates for plot
start_date <- as.Date("2001-01-01")
end_date <- as.Date("2015-01-01")
# Create sequence of dates for x-axis
dates_seq <- seq(from = start_date, to = end_date, by = "year")
# Subset event_df to only include dates after or equal to start_date
event_df_sub <- event_df[dates >= start_date, ]
pl <- ggplot(event_df, aes(x = dates, y = coef)) +
geom_point() +geom_line()+
#geom_errorbar(aes(ymin = ci_lower, ymax = ci_upper), width = 0.2) +
geom_ribbon(aes(ymin = ci_lower, ymax = ci_upper), alpha = 0.2) +
geom_hline(yintercept = 0, linetype = "dashed", color = "red") +
geom_vline(xintercept = as.Date("2005-01-01"), linetype = "dashed", color = "blue") +
scale_y_continuous(labels = function(x) paste0(x * 100, "%")) +
scale_x_date(date_labels = "%Y", limits = c(start_date, end_date), breaks = "5 year") +
labs(title = "Panel A: Rent 0", y = "Percentage")+ #title = "Fixed effect of RENT_0") +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
plot.title = element_text(size = 14, face = "bold"),
axis.title.x = element_blank(),
axis.title.y = element_text(size = 10, face = "bold"),
axis.text = element_text(size = 12),
legend.text = element_text(size = 12),
plot.caption = element_text(size = 12, face = c("italic", "bold")),
axis.line = element_line(color = "black", size = 0.5),
panel.background = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_line(color = "black", size = 0.1))## Warning: The `size` argument of `element_line()` is deprecated as of ggplot2 3.4.0.
## ℹ Please use the `linewidth` argument instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.## Warning: Vectorized input to `element_text()` is not officially supported.
## ℹ Results may be unexpected or may change in future versions of ggplot2.pl
Rent 1
cse.01 <- conleyreg(log(RENT_1) ~ Price + Unemployee.Rate +
Personal.Income + Resident.Population +
treat:(d.2001 + d.2002 + d.2003 +
d.2005 + d.2006 + d.2007 + d.2008 +
d.2009 + d.2010 + d.2011 + d.2012 +
d.2013 + d.2014 + d.2015) | MSA + YEAR,
data = df,
dist_cutoff = 300,
st_distance = TRUE,
lat = "Y",
lon = "X",
)## Checking spatial attributes## Estimating model## Warning: In fixest_env(fml = fml, data = data, weights = weig...:
## OpenMP not detected: cannot use 4 threads, single-threaded mode instead.## Estimating distance matrix and addressing spatial correlationsummary(cse.01)## Estimate Std. Error t value Pr(>|t|)
## Min. :-0.001283 Min. :3.260e-06 Min. :-2.106 Min. :0.0000000
## 1st Qu.: 0.007848 1st Qu.:5.125e-03 1st Qu.: 2.080 1st Qu.:0.0002195
## Median : 0.107106 Median :2.970e-02 Median : 3.077 Median :0.0024954
## Mean : 0.080398 Mean :2.420e-02 Mean : 2.813 Mean :0.0699408
## 3rd Qu.: 0.134174 3rd Qu.:3.881e-02 3rd Qu.: 3.708 3rd Qu.:0.0291474
## Max. : 0.173091 Max. :4.870e-02 Max. : 5.580 Max. :0.9504810colnames(cse.01)[2] <- "Conley SE"
coef_df1 <- cse.01[,] %>%
as_tibble() %>%
mutate(variables = rownames(cse.01))
coef_df1 <- coef_df1[5:18,]
event_df1 <- data.frame(time = 1:nrow(coef_df1),
coef = coef_df1$Estimate,
cse = coef_df1$`Conley SE`)
event_df1 <- rbind(event_df1, c("time" = 3.5, "coef" = 0, "cse" = 0))
event_df1 <- arrange(event_df1, time)
# Calculate confidence intervals
event_df1$ci_upper <- event_df1$coef + 1.96 * event_df1$cse
event_df1$ci_lower <- event_df1$coef - 1.96 * event_df1$cse
# Subset event_df to only include dates after or equal to start_date
event_df_sub1 <- event_df1[dates >= start_date, ]
pl1 <- ggplot(event_df1, aes(x = dates, y = coef)) +
geom_point() +geom_line()+
#geom_errorbar(aes(ymin = ci_lower, ymax = ci_upper), width = 0.2) +
geom_ribbon(aes(ymin = ci_lower, ymax = ci_upper), alpha = 0.2) +
geom_hline(yintercept = 0, linetype = "dashed", color = "red") +
geom_vline(xintercept = as.Date("2005-01-01"), linetype = "dashed", color = "blue") +
scale_y_continuous(labels = function(x) paste0(x * 100, "%")) +
scale_x_date(date_labels = "%Y", limits = c(start_date, end_date), breaks = "5 year") +
labs(title = "Panel B: Rent 1", y = "Percentage")+ #title = "Fixed effect of RENT_1") +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
plot.title = element_text(size = 14, face = "bold"),
axis.title.x = element_blank(),
axis.title.y = element_text(size = 10, face = "bold"),
axis.text = element_text(size = 12),
legend.text = element_text(size = 12),
plot.caption = element_text(size = 12, face = c("italic", "bold")),
axis.line = element_line(color = "black", size = 0.5),
panel.background = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_line(color = "black", size = 0.1))## Warning: Vectorized input to `element_text()` is not officially supported.
## ℹ Results may be unexpected or may change in future versions of ggplot2.pl1
Rent 2
cse.02 <- conleyreg(log(RENT_2) ~ Price + Unemployee.Rate +
Personal.Income + Resident.Population +
treat:(d.2001 + d.2002 + d.2003 +
d.2005 + d.2006 + d.2007 + d.2008 +
d.2009 + d.2010 + d.2011 + d.2012 +
d.2013 + d.2014 + d.2015) | MSA + YEAR,
data = df,
dist_cutoff = 300,
st_distance = TRUE,
lat = "Y",
lon = "X",
)## Checking spatial attributes## Estimating model## Warning: In fixest_env(fml = fml, data = data, weights = weig...:
## OpenMP not detected: cannot use 4 threads, single-threaded mode instead.## Estimating distance matrix and addressing spatial correlationsummary(cse.02)## Estimate Std. Error t value Pr(>|t|)
## Min. :-0.0001241 Min. :2.910e-06 Min. :-2.733 Min. :0.0000004
## 1st Qu.: 0.0082509 1st Qu.:4.928e-03 1st Qu.: 2.228 1st Qu.:0.0011086
## Median : 0.0961914 Median :2.992e-02 Median : 2.511 Median :0.0089079
## Mean : 0.0704461 Mean :2.420e-02 Mean : 2.452 Mean :0.0800094
## 3rd Qu.: 0.1127988 3rd Qu.:3.861e-02 3rd Qu.: 3.273 3rd Qu.:0.0236205
## Max. : 0.1582695 Max. :4.789e-02 Max. : 5.091 Max. :0.8778430colnames(cse.02)[2] <- "Conley SE"
coef_df2 <- cse.02[,] %>%
as_tibble() %>%
mutate(variables = rownames(cse.02))
coef_df2 <- coef_df2[5:18,]
event_df2 <- data.frame(time = 1:nrow(coef_df2),
coef = coef_df2$Estimate,
cse = coef_df2$`Conley SE`)
event_df2 <- rbind(event_df2, c("time" = 3.5, "coef" = 0, "cse" = 0))
event_df2 <- arrange(event_df2, time)
# Calculate confidence intervals
event_df2$ci_upper <- event_df2$coef + 1.96 * event_df2$cse
event_df2$ci_lower <- event_df2$coef - 1.96 * event_df2$cse
# Subset event_df to only include dates after or equal to start_date
event_df_sub2 <- event_df2[dates >= start_date, ]
pl2 <- ggplot(event_df2, aes(x = dates, y = coef)) +
geom_point() +
geom_line() +
#geom_errorbar(aes(ymin = ci_lower, ymax = ci_upper), width = 0.2) +
geom_ribbon(aes(ymin = ci_lower, ymax = ci_upper), alpha = 0.2) +
geom_hline(yintercept = 0, linetype = "dashed", color = "red") +
geom_vline(xintercept = as.Date("2005-01-01"), linetype = "dashed", color = "blue") +
scale_y_continuous(labels = function(x) paste0(x * 100, "%")) +
scale_x_date(date_labels = "%Y", limits = c(start_date, end_date), breaks = "5 year") +
labs(title = "Panel C: RENT 2", y = "Percentage")+ #title = "Fixed effect of RENT_1") +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
plot.title = element_text(size = 14, face = "bold"),
axis.title.x = element_blank(),
axis.title.y = element_text(size = 10, face = "bold"),
axis.text = element_text(size = 12),
legend.text = element_text(size = 12),
plot.caption = element_text(size = 12, face = c("italic", "bold")),
axis.line = element_line(color = "black", size = 0.5),
panel.background = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_line(color = "black", size = 0.1))## Warning: Vectorized input to `element_text()` is not officially supported.
## ℹ Results may be unexpected or may change in future versions of ggplot2.pl2
Rent 3
cse.03 <- conleyreg(log(RENT_3) ~ Price + Unemployee.Rate +
Personal.Income + Resident.Population +
treat:(d.2001 + d.2002 + d.2003 +
d.2005 + d.2006 + d.2007 + d.2008 +
d.2009 + d.2010 + d.2011 + d.2012 +
d.2013 + d.2014 + d.2015) | MSA + YEAR,
data = df,
dist_cutoff = 300,
st_distance = TRUE,
lat = "Y",
lon = "X",
)## Checking spatial attributes## Estimating model## Warning: In fixest_env(fml = fml, data = data, weights = weig...:
## OpenMP not detected: cannot use 4 threads, single-threaded mode instead.## Estimating distance matrix and addressing spatial correlationsummary(cse.03)## Estimate Std. Error t value Pr(>|t|)
## Min. :-0.048740 Min. :3.280e-06 Min. :-3.7820 Min. :0.0001095
## 1st Qu.: 0.006022 1st Qu.:5.236e-03 1st Qu.: 0.8093 1st Qu.:0.0098364
## Median : 0.032740 Median :3.101e-02 Median : 1.4577 Median :0.0647285
## Mean : 0.035300 Mean :2.483e-02 Mean : 1.2705 Mean :0.1689422
## 3rd Qu.: 0.057885 3rd Qu.:4.031e-02 3rd Qu.: 2.5829 3rd Qu.:0.2682266
## Max. : 0.105041 Max. :5.018e-02 Max. : 3.8827 Max. :0.7951908colnames(cse.03)[2] <- "Conley SE"
coef_df3 <- cse.03[,] %>%
as_tibble() %>%
mutate(variables = rownames(cse.03))
coef_df3 <- coef_df3[5:18,]
event_df3 <- data.frame(time = 1:nrow(coef_df3),
coef = coef_df3$Estimate,
cse = coef_df3$`Conley SE`)
event_df3 <- rbind(event_df3, c("time" = 3.5, "coef" = 0, "cse" = 0))
event_df3 <- arrange(event_df3, time)
# Calculate confidence intervals
event_df3$ci_upper <- event_df3$coef + 1.96 * event_df3$cse
event_df3$ci_lower <- event_df3$coef - 1.96 * event_df3$cse
# Subset event_df to only include dates after or equal to start_date
event_df_sub3 <- event_df3[dates >= start_date, ]
pl3 <- ggplot(event_df3, aes(x = dates, y = coef)) +
geom_point() +
geom_line() +
#geom_errorbar(aes(ymin = ci_lower, ymax = ci_upper), width = 0.2) +
geom_ribbon(aes(ymin = ci_lower, ymax = ci_upper), alpha = 0.2) +
geom_hline(yintercept = 0, linetype = "dashed", color = "red") +
geom_vline(xintercept = as.Date("2005-01-01"), linetype = "dashed", color = "blue") +
scale_y_continuous(labels = function(x) paste0(x * 100, "%")) +
scale_x_date(date_labels = "%Y", limits = c(start_date, end_date), breaks = "5 year") +
labs(title = "Panel D: RENT 3", y = "Percentage")+ #title = "Fixed effect of RENT_1") +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
plot.title = element_text(size = 14, face = "bold"),
axis.title.x = element_blank(),
axis.title.y = element_text(size = 10, face = "bold"),
axis.text = element_text(size = 12),
legend.text = element_text(size = 12),
plot.caption = element_text(size = 12, face = c("italic", "bold")),
axis.line = element_line(color = "black", size = 0.5),
panel.background = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_line(color = "black", size = 0.1))## Warning: Vectorized input to `element_text()` is not officially supported.
## ℹ Results may be unexpected or may change in future versions of ggplot2.pl3
Rent 4
cse.04 <- conleyreg(log(RENT_4) ~ Price + Unemployee.Rate +
Personal.Income + Resident.Population +
treat:(d.2001 + d.2002 + d.2003 +
d.2005 + d.2006 + d.2007 + d.2008 +
d.2009 + d.2010 + d.2011 + d.2012 +
d.2013 + d.2014 + d.2015) | MSA + YEAR,
data = df,
dist_cutoff = 300,
st_distance = TRUE,
lat = "Y",
lon = "X",
)## Checking spatial attributes## Estimating model## Warning: In fixest_env(fml = fml, data = data, weights = weig...:
## OpenMP not detected: cannot use 4 threads, single-threaded mode instead.## Estimating distance matrix and addressing spatial correlationsummary(cse.04)## Estimate Std. Error t value Pr(>|t|)
## Min. :-0.031638 Min. :2.820e-06 Min. :-1.8679 Min. :0.000106
## 1st Qu.: 0.001889 1st Qu.:5.830e-03 1st Qu.: 0.8029 1st Qu.:0.056766
## Median : 0.035738 Median :4.352e-02 Median : 1.0744 Median :0.199581
## Mean : 0.036327 Mean :3.319e-02 Mean : 1.2492 Mean :0.215039
## 3rd Qu.: 0.061519 3rd Qu.:5.294e-02 3rd Qu.: 1.8707 3rd Qu.:0.336023
## Max. : 0.100210 Max. :5.833e-02 Max. : 3.8909 Max. :0.488579colnames(cse.04)[2] <- "Conley SE"
coef_df4 <- cse.04[,] %>%
as_tibble() %>%
mutate(variables = rownames(cse.04))
coef_df4 <- coef_df4[5:18,]
event_df4 <- data.frame(time = 1:nrow(coef_df4),
coef = coef_df4$Estimate,
cse = coef_df4$`Conley SE`)
event_df4 <- rbind(event_df4, c("time" = 3.5, "coef" = 0, "cse" = 0))
event_df4 <- arrange(event_df4, time)
# Calculate confidence intervals
event_df4$ci_upper <- event_df4$coef + 1.96 * event_df4$cse
event_df4$ci_lower <- event_df4$coef - 1.96 * event_df4$cse
# Subset event_df to only include dates after or equal to start_date
event_df_sub4 <- event_df4[dates >= start_date, ]
pl4 <- ggplot(event_df4, aes(x = dates, y = coef)) +
geom_point() +
geom_line() +
#geom_errorbar(aes(ymin = ci_lower, ymax = ci_upper), width = 0.2) +
geom_ribbon(aes(ymin = ci_lower, ymax = ci_upper), alpha = 0.2) +
geom_hline(yintercept = 0, linetype = "dashed", color = "red") +
geom_vline(xintercept = as.Date("2005-01-01"), linetype = "dashed", color = "blue") +
scale_y_continuous(labels = function(x) paste0(x * 100, "%")) +
scale_x_date(date_labels = "%Y", limits = c(start_date, end_date), breaks = "5 year") +
labs(title = "Panel E: RENT 4", y = "Percentage") + #title = "Fixed effect of RENT_1") +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
plot.title = element_text(size = 14, face = "bold"),
axis.title.x = element_blank(),
axis.title.y = element_text(size = 10, face = "bold"),
axis.text = element_text(size = 12),
legend.text = element_text(size = 12),
plot.caption = element_text(size = 12, face = c("italic", "bold")),
axis.line = element_line(color = "black", size = 0.5),
panel.background = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_line(color = "black", size = 0.1))## Warning: Vectorized input to `element_text()` is not officially supported.
## ℹ Results may be unexpected or may change in future versions of ggplot2.pl4
House Price Index
cse.prc <- conleyreg(log(Price) ~ Unemployee.Rate +
Personal.Income + Resident.Population +
treat:(d.2001 + d.2002 + d.2003 +
d.2005 + d.2006 + d.2007 + d.2008 +
d.2009 + d.2010 + d.2011 + d.2012 +
d.2013 + d.2014 + d.2015) | MSA + YEAR,
data = df,
dist_cutoff = 300,
st_distance = TRUE,
lat = "Y",
lon = "X",
)## Checking spatial attributes## Estimating model## Warning: In fixest_env(fml = fml, data = data, weights = weig...:
## OpenMP not detected: cannot use 4 threads, single-threaded mode instead.## Estimating distance matrix and addressing spatial correlationsummary(cse.prc)## Estimate Std. Error t value Pr(>|t|)
## Min. :-0.10543 Min. :4.510e-06 Min. :-3.915 Min. :0.0000001
## 1st Qu.:-0.02722 1st Qu.:2.051e-02 1st Qu.:-1.068 1st Qu.:0.0000678
## Median : 0.02609 Median :2.548e-02 Median : 2.217 Median :0.0028688
## Mean : 0.03769 Mean :2.323e-02 Mean : 1.362 Mean :0.1187125
## 3rd Qu.: 0.10499 3rd Qu.:3.019e-02 3rd Qu.: 3.999 3rd Qu.:0.2167757
## Max. : 0.19225 Max. :3.550e-02 Max. : 5.415 Max. :0.9059312colnames(cse.prc)[2] <- "Conley SE"
coef_df5 <- cse.prc[,] %>%
as_tibble() %>%
mutate(variables = rownames(cse.prc))
coef_df5 <- coef_df5[5:18,]
event_df5 <- data.frame(time = 1:nrow(coef_df5),
coef = coef_df5$Estimate,
cse = coef_df5$`Conley SE`)
event_df5 <- rbind(event_df5, c("time" = 3.5, "coef" = 0, "cse" = 0))
event_df5 <- arrange(event_df5, time)
# Calculate confidence intervals
event_df5$ci_upper <- event_df5$coef + 1.96 * event_df5$cse
event_df5$ci_lower <- event_df5$coef - 1.96 * event_df5$cse
# Subset event_df to only include dates after or equal to start_date
event_df_sub5 <- event_df5[dates >= start_date, ]
pl5 <- ggplot(event_df5, aes(x = dates, y = coef)) +
geom_point() +
geom_line() +
#geom_errorbar(aes(ymin = ci_lower, ymax = ci_upper), width = 0.2) +
geom_ribbon(aes(ymin = ci_lower, ymax = ci_upper), alpha = 0.2) +
geom_hline(yintercept = 0, linetype = "dashed", color = "red") +
geom_vline(xintercept = as.Date("2005-01-01"), linetype = "dashed", color = "blue") +
scale_y_continuous(labels = function(x) paste0(x * 100, "%")) +
scale_x_date(date_labels = "%Y", limits = c(start_date, end_date), breaks = "5 year") +
labs(title = "Panel F: House Price Index", y = "Percentage") + #title = "Fixed effect of RENT_1") +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
plot.title = element_text(size = 14, face = "bold"),
axis.title.x = element_blank(),
axis.title.y = element_text(size = 10, face = "bold"),
axis.text = element_text(size = 12),
legend.text = element_text(size = 12),
plot.caption = element_text(size = 12, face = c("italic", "bold")),
axis.line = element_line(color = "black", size = 0.5),
panel.background = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_line(color = "black", size = 0.1))## Warning: Vectorized input to `element_text()` is not officially supported.
## ℹ Results may be unexpected or may change in future versions of ggplot2.pl5## Warning: Removed 1 rows containing missing values (`geom_point()`).## Warning: Removed 1 row containing missing values (`geom_line()`).
All Plots
grid.arrange(pl, pl1, pl2, pl3, pl4, pl5, nrow = 2)## Warning: Removed 1 rows containing missing values (`geom_point()`).## Warning: Removed 1 row containing missing values (`geom_line()`).