assignment_3

Author

Sam Kuhn

Part 1: Basic Fixed Effects model

## Load data
#| echo: false
#| message: false
library(tidyverse)
library(haven)
library(gt)
library(modelsummary)
library(fixest)
library(plm)
library(kableExtra)

rent_df <- haven::read_dta(here::here("data", "rent.dta"))
davis_df <- haven::read_dta(here::here("data", "davis2004.dta"))

1. Pooled OLS model

Pooled OLS

pooled_model <- plm(lnrent ~ y90 + lnpopulation + lninc + studentpop,
  model = "pooling",
  data = rent_df
)
kable(broom::tidy(pooled_model), digits = 3, caption = "Pooled Model")

Pooled Model
term	estimate	std.error	statistic	p.value
(Intercept)	-0.569	0.535	-1.063	0.290
y90	0.262	0.035	7.543	0.000
lnpopulation	0.041	0.023	1.807	0.073
lninc	0.571	0.053	10.762	0.000
studentpop	0.005	0.001	4.949	0.000

4. Estimate demeaned FE model

Demeaned model

#calculate city-level averages for all variables, and subtract from observed values
rent_df <- rent_df |>
  group_by(city) |>
  mutate(across(!matches(c("year", "y90")), ~ .x - mean(.x), .names = "demeaned_{col}"))

#If year == 80 for demeaned vars, replace with NA
rent_df <- rent_df |> 
  mutate(across(starts_with("demeaned"), ~ replace(.x, year == 80, NA)))

demeaned_model <- lm(demeaned_rent ~ y90 + demeaned_population + demeaned_studentpop, data = rent_df)
summary(demeaned_model)


Call:
lm(formula = demeaned_rent ~ y90 + demeaned_population + demeaned_studentpop, 
    data = rent_df)

Residuals:
    Min      1Q  Median      3Q     Max 
-37.870 -20.293  -7.794  12.492  99.957 

Coefficients: (1 not defined because of singularities)
                     Estimate Std. Error t value Pr(>|t|)    
(Intercept)         8.326e+01  4.678e+00  17.796   <2e-16 ***
y90                        NA         NA      NA       NA    
demeaned_population 1.350e-03  5.251e-04   2.571   0.0126 *  
demeaned_studentpop 5.536e+00  2.763e+00   2.004   0.0496 *  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 30.53 on 61 degrees of freedom
  (64 observations deleted due to missingness)
Multiple R-squared:  0.1428,    Adjusted R-squared:  0.1147 
F-statistic: 5.083 on 2 and 61 DF,  p-value: 0.009087

6. Estimate FE model to verify identical estimates to (4)

demeaned_model_check <- plm(lnrent ~ y90 + lnpopulation + lninc + studentpop,
  model = "within",
  index = c("city", "year"),
  data = rent_df
)
kable(broom::tidy(demeaned_model_check), digits = 3, caption = "Pooled Model")

Pooled Model
term	estimate	std.error	statistic	p.value
y90	0.386	0.037	10.469	0.000
lnpopulation	0.072	0.088	0.818	0.417
lninc	0.310	0.066	4.663	0.000
studentpop	0.011	0.004	2.711	0.009

Part 2: Advanced FE model

1. What are general issues apply housing values for estimating the value of health risks?

2. Run the following model with using the samples from both Churchill and Lyon County.

\(ln(price) = \alpha + \beta CancerRisk+ \epsilon\)

churchill_lyon_model <- fixest::feols(
  lnprice ~ cancerRisk,
  data = davis_df
)
modelsummary(churchill_lyon_model, vcov = cluster ~ group, stars = TRUE)

	(1)
(Intercept)	11.625***
	(0.006)
cancerRisk	−0.108***
	(0.019)
Num.Obs.	10204
R2	0.003
R2 Adj.	0.003
AIC	9887.9
BIC	9902.3
RMSE	0.39
Std.Errors	by: cluster
+ p < 0.1, * p < 0.05, p < 0.01, * p < 0.001

3. Run the following model with using the sample only from Churchill County.

\(ln(price) = \alpha + \beta CancerRisk+ \epsilon\)

churchill_model <- davis_df |>
  filter(churchill == 1) |>
  fixest::feols(lnprice ~ cancerRisk)
modelsummary(churchill_model, vcov = cluster ~ group, stars = TRUE)

	(1)
(Intercept)	11.588***
	(0.009)
cancerRisk	−0.064**
	(0.020)
Num.Obs.	3596
R2	0.002
R2 Adj.	0.002
AIC	3754.5
BIC	3766.9
RMSE	0.41
Std.Errors	by: cluster
+ p < 0.1, * p < 0.05, p < 0.01, * p < 0.001

4. Run the following model with using the samples from both Churchill and Lyon County.

\(ln(price) = \alpha + \beta_1 CancerRisk + \beta_2 Churchill + \epsilon\)

churchill_lyon_fe_model <- fixest::feols(
  lnprice ~ cancerRisk | churchill,
  data = davis_df
)
modelsummary(churchill_lyon_fe_model, vcov = cluster ~ group, stars = TRUE)

	(1)
cancerRisk	−0.064
	(0.020)
Num.Obs.	10204
R2	0.007
R2 Adj.	0.006
R2 Within	0.001
R2 Within Adj.	0.001
AIC	9852.8
BIC	9874.5
RMSE	0.39
Std.Errors	by: cluster
FE: churchill	X
+ p < 0.1, * p < 0.05, p < 0.01, * p < 0.001

5. Run the following model with using the samples from both Churchill and Lyon County.

\(ln(price) = \alpha + \beta_1 CancerRisk + \beta_2 Churchill + \beta_3X + \mu_y + \gamma_m + \epsilon\)

full_fe_model <- fixest::feols(
  lnprice ~ cancerRisk + lotsize + lotsizesq + floorsize + age + agesq + cl_1 + cl_2 + cl_3 + cl_4 + cl_5 + cl_6 + cl_7 + cl_8 + cl_9 | churchill + year + month,
  data = davis_df
)

The variable 'cl_9' has been removed because of collinearity (see $collin.var).

modelsummary(full_fe_model, vcov = cluster ~ group, stars = TRUE)

	(1)
cancerRisk	−0.156+
	(0.017)
lotsize	0.012+
	(0.002)
lotsizesq	0.000
	(0.000)
floorsize	0.044*
	(0.001)
age	−0.006+
	(0.001)
agesq	0.000
	(0.000)
cl_1	−0.480+
	(0.053)
cl_2	−0.606+
	(0.055)
cl_3	−0.621+
	(0.056)
cl_4	−0.499+
	(0.050)
cl_5	−0.405+
	(0.049)
cl_6	−0.288
	(0.048)
cl_7	−0.223
	(0.046)
cl_8	−0.120
	(0.060)
Num.Obs.	10204
R2	0.633
R2 Adj.	0.632
R2 Within	0.624
R2 Within Adj.	0.624
AIC	−244.8
BIC	37.2
RMSE	0.24
Std.Errors	by: cluster
FE: churchill	X
FE: year	X
FE: month	X
+ p < 0.1, * p < 0.05, p < 0.01, * p < 0.001

6. Run the following model with using the samples from both Churchill and Lyon County and with houses that were sold multiple times in the dataset.

\(ln(price) = \alpha + \beta_1 CancerRisk + \mu_y + \gamma_m + \pi_h + \epsilon\)

davis_df |> 
  group_by(houseID) |> 
  summarise(count_sale = n()) |> 
  arrange(desc(count_sale))

# A tibble: 7,450 × 2
     houseID count_sale
       <dbl>      <int>
 1 300121502          7
 2 300138402          6
 3 200120118          5
 4 200145201          5
 5 200163112          5
 6 200725135          5
 7 200875321          5
 8 201006115          5
 9 301946402          5
10 301966308          5
# ℹ 7,440 more rows

davis_df |> 
  group_by(houseID) |> 
  summarise(count_sale = n()) |> 
  ggplot(aes(count_sale)) +
  geom_histogram(bins = 10) +
  labs(
    title = "Distribution of houses sold more than once"
  ) +
  theme_minimal()

house_fe_model <- fixest::feols(
  lnprice ~ cancerRisk | year + month + houseID,
  data = davis_df
)
modelsummary(house_fe_model, vcov = cluster ~ group, stars = TRUE)

	(1)
cancerRisk	−0.140***
	(0.022)
Num.Obs.	10204
R2	0.944
R2 Adj.	0.791
R2 Within	0.017
R2 Within Adj.	0.017
AIC	−4586.0
BIC	49455.0
RMSE	0.09
Std.Errors	by: cluster
FE: year	X
FE: month	X
FE: houseID	X
+ p < 0.1, * p < 0.05, p < 0.01, * p < 0.001