Regression Discontinuity in Time Analysis: NQI Impact on Patenting

Setup and Data Loading

# Load required packages
library(tidyverse)
library(lubridate)
library(fixest)
library(modelsummary)
library(marginaleffects)
library(kableExtra)

# Set theme for all plots
theme_set(theme_minimal(base_size = 12))

# Load data
df <- read_csv("/Users/amitdas/Downloads/UO_phd/data_constructed/QC_Emergence/data/data_constructed/total_qc_patents_17Nov25.csv")

# Parse dates and create time variables
df <- df %>%
  mutate(
    patent_date = ymd(patent_date),
    patent_year = year(patent_date),
    patent_qtr = quarter(patent_date),
    yq = paste0(patent_year, "Q", patent_qtr),
    lobbied_status = ifelse(qcLobbied == 1, 1, 0)
  )

# Data summary
cat("=== DATA SUMMARY ===\n")

=== DATA SUMMARY ===

cat("Unique patents:", n_distinct(df$patent_id), "\n")

Unique patents: 5412 

cat("Unique firms:", n_distinct(df$cleaned_Firms_patented), "\n")

Unique firms: 655 

cat("Date range:", as.character(min(df$patent_date)), "to", 
    as.character(max(df$patent_date)), "\n")

Date range: 2000-01-04 to 2024-12-31 

cat("Patents from lobbied firms:", sum(df$lobbied_status == 1), "\n")

Patents from lobbied firms: 2500 

cat("Patents from non-lobbied firms:", sum(df$lobbied_status == 0), "\n")

Patents from non-lobbied firms: 3044 

Descriptive Analysis

Figure 1: Quarterly Patent Counts Over Time

# Aggregate to quarterly level
qtr_counts <- df %>%
  group_by(patent_year, patent_qtr, yq) %>%
  summarise(unique_patents = n_distinct(patent_id), .groups = "drop") %>%
  mutate(
    time = patent_year + (patent_qtr - 1) / 4,
    period = case_when(
      patent_year <= 2017 ~ "before",
      patent_year == 2018 ~ "2018",
      patent_year >= 2019 ~ "after"
    )
  )

# Plot excluding 2018
qtr_counts %>%
  filter(patent_year != 2018) %>%
  ggplot(aes(x = time, y = unique_patents)) +
  geom_point(aes(color = period), size = 3, alpha = 0.7) +
  geom_smooth(
    data = . %>% filter(period == "before"),
    method = "lm", se = TRUE,
    color = "#E41A1C", fill = "#E41A1C", alpha = 0.2, linewidth = 1.2
  ) +
  geom_smooth(
    data = . %>% filter(period == "after"),
    method = "lm", se = TRUE,
    color = "#377EB8", fill = "#377EB8", alpha = 0.2, linewidth = 1.2
  ) +
  geom_vline(xintercept = 2018, linetype = "dashed", color = "gray50", linewidth = 0.8) +
  annotate("text", x = 2018, y = max(qtr_counts$unique_patents[qtr_counts$patent_year != 2018]) * 0.95, 
           label = "2018\n(NQI)", size = 3.5, color = "gray30") +
  scale_color_manual(
    values = c("before" = "#E41A1C", "after" = "#377EB8"),
    labels = c("before" = "2000-2017", "after" = "2019-2024")
  ) +
  labs(
    title = "Quarterly Patent Counts Over Time",
    subtitle = "Unique patent counts with trend lines before and after 2018 NQI policy",
    x = "Year",
    y = "Number of Unique Patents",
    color = "Period"
  ) +
  theme(
    plot.title = element_text(face = "bold", size = 14),
    plot.subtitle = element_text(size = 11, color = "gray40"),
    legend.position = "bottom",
    panel.grid.minor = element_blank()
  )

Figure 1

Figure 2: Patent Counts by Lobbying Status

# Aggregate by quarter and lobbying status
qtr_counts_lobby <- df %>%
  group_by(patent_year, patent_qtr, yq, lobbied_status) %>%
  summarise(unique_patents = n_distinct(patent_id), .groups = "drop") %>%
  mutate(
    time = patent_year + (patent_qtr - 1) / 4,
    period = case_when(
      patent_year <= 2017 ~ "before",
      patent_year == 2018 ~ "2018",
      patent_year >= 2019 ~ "after"
    ),
    lobby_label = ifelse(lobbied_status == 1, "Lobbied", "Not Lobbied")
  )

# Plot with facets
qtr_counts_lobby %>%
  filter(patent_year != 2018) %>%
  ggplot(aes(x = time, y = unique_patents)) +
  geom_point(aes(color = period), size = 3, alpha = 0.7) +
  geom_smooth(
    data = . %>% filter(period == "before"),
    method = "lm", se = TRUE,
    color = "#E41A1C", fill = "#E41A1C", alpha = 0.2, linewidth = 1.2
  ) +
  geom_smooth(
    data = . %>% filter(period == "after"),
    method = "lm", se = TRUE,
    color = "#377EB8", fill = "#377EB8", alpha = 0.2, linewidth = 1.2
  ) +
  geom_vline(xintercept = 2018, linetype = "dashed", color = "gray50", linewidth = 0.8) +
  facet_wrap(~lobby_label, ncol = 2) +
  scale_color_manual(
    values = c("before" = "#E41A1C", "after" = "#377EB8"),
    labels = c("before" = "2000-2017", "after" = "2019-2024")
  ) +
  labs(
    title = "Quarterly Patent Counts: Lobbied vs Non-Lobbied Firms",
    subtitle = "Trend lines before and after 2018 NQI policy",
    x = "Year",
    y = "Number of Unique Patents",
    color = "Period"
  ) +
  theme(
    plot.title = element_text(face = "bold", size = 14),
    plot.subtitle = element_text(size = 11, color = "gray40"),
    legend.position = "bottom",
    panel.grid.minor = element_blank(),
    strip.text = element_text(face = "bold", size = 12)
  )

Figure 2

# Calculate and display trend slopes
cat("\n=== TREND LINE SLOPES ===\n")

=== TREND LINE SLOPES ===

plot_data <- qtr_counts_lobby %>% filter(patent_year != 2018)

# Before 2018
before_data <- plot_data %>% filter(period == "before")
lm_before_lobby <- lm(unique_patents ~ time, data = filter(before_data, lobbied_status == 1))
lm_before_notlobby <- lm(unique_patents ~ time, data = filter(before_data, lobbied_status == 0))

cat("Before 2018 (2000-2017):\n")

Before 2018 (2000-2017):

cat("  Lobbied slope:", round(coef(lm_before_lobby)[2], 3), "\n")

  Lobbied slope: 0.761 

cat("  Not Lobbied slope:", round(coef(lm_before_notlobby)[2], 3), "\n\n")

  Not Lobbied slope: 0.799 

# After 2018
after_data <- plot_data %>% filter(period == "after")
lm_after_lobby <- lm(unique_patents ~ time, data = filter(after_data, lobbied_status == 1))
lm_after_notlobby <- lm(unique_patents ~ time, data = filter(after_data, lobbied_status == 0))

cat("After 2018 (2019-2024):\n")

After 2018 (2019-2024):

cat("  Lobbied slope:", round(coef(lm_after_lobby)[2], 3), "\n")

  Lobbied slope: 14.981 

cat("  Not Lobbied slope:", round(coef(lm_after_notlobby)[2], 3), "\n")

  Not Lobbied slope: 9.263 

Regression Discontinuity in Time (RDiT) Analysis

Data Preparation

# Create firm-quarter panel
firm_qtr_panel <- df %>%
  group_by(firm = cleaned_Firms_patented, patent_year, patent_qtr, lobbied_status) %>%
  summarise(patent_count = n_distinct(patent_id), .groups = "drop") %>%
  mutate(
    yq = paste0(patent_year, "Q", patent_qtr),
    # Continuous time variable (quarters since 2000Q1)
    time = (patent_year - 2000) * 4 + patent_qtr,
    # Time centered at 2018Q1 (cutoff = 73)
    cutoff = (2018 - 2000) * 4 + 1,
    Duration = time - cutoff,
    # Post-policy dummy
    PostPolicy = ifelse(patent_year >= 2019, 1, 0)
  )

# Exclude 2018 for cleaner RDiT
rdit_data <- firm_qtr_panel %>%
  filter(patent_year != 2018) %>%
  mutate(
    Duration1 = Duration,
    Duration2 = Duration^2,
    Duration3 = Duration^3
  )

cat("=== PANEL DATA SUMMARY ===\n")

=== PANEL DATA SUMMARY ===

cat("Total firm-quarter observations:", nrow(rdit_data), "\n")

Total firm-quarter observations: 2882 

cat("Unique firms:", n_distinct(rdit_data$firm), "\n")

Unique firms: 643 

cat("Time range:", min(rdit_data$patent_year), "to", max(rdit_data$patent_year), "\n")

Time range: 2000 to 2024 

Optimal Polynomial Order Selection

# Test polynomial orders 1-4
model_fit <- tibble(poly_order = 1:4, AIC = NA_real_, BIC = NA_real_)

for (p in 1:4) {
  temp_data <- rdit_data
  for (i in 1:p) {
    temp_data[[paste0("Duration", i)]] <- temp_data$Duration^i
  }
  
  formula_str <- paste0("patent_count ~ PostPolicy + ",
                       paste0("Duration", 1:p, collapse = " + "))
  model <- lm(as.formula(formula_str), data = temp_data)
  
  model_fit$AIC[p] <- AIC(model)
  model_fit$BIC[p] <- BIC(model)
}

cat("=== MODEL FIT COMPARISON ===\n")

=== MODEL FIT COMPARISON ===

print(model_fit)

# A tibble: 4 × 3
  poly_order    AIC    BIC
       <int>  <dbl>  <dbl>
1          1 13599. 13623.
2          2 13596. 13626.
3          3 13598. 13634.
4          4 13599. 13641.

cat("\nOptimal polynomial order (by AIC):", which.min(model_fit$AIC), "\n")

Optimal polynomial order (by AIC): 2 

cat("Optimal polynomial order (by BIC):", which.min(model_fit$BIC), "\n")

Optimal polynomial order (by BIC): 1 

cat("\nUsing polynomial order 3 for analysis\n")

Using polynomial order 3 for analysis

Hypothesis 1: Main Effect of NQI 2018

H1: NQI 2018 will increase patenting counts of firms

# Model 1a: Simple RDiT (no fixed effects)
h1_simple <- lm(
  patent_count ~ PostPolicy + Duration1 + Duration2 + Duration3,
  data = rdit_data
)

# Model 1b: RDiT with firm fixed effects (PREFERRED)
h1_fe <- feols(
  patent_count ~ PostPolicy + Duration1 + Duration2 + Duration3 | firm,
  data = rdit_data,
  cluster = ~ firm
)

# Model 1c: RDiT with different slopes before/after
h1_interact <- feols(
  patent_count ~ PostPolicy + Duration1 + Duration2 + Duration3 +
    PostPolicy:Duration1 + PostPolicy:Duration2 + PostPolicy:Duration3 | firm,
  data = rdit_data,
  cluster = ~ firm
)

Table 1: H1 Regression Results

# Create publication-quality table
modelsummary(
  list(
    "Simple RDiT" = h1_simple,
    "Firm FE" = h1_fe,
    "Different Slopes" = h1_interact
  ),
  stars = c('*' = .1, '**' = .05, '***' = .01),
  coef_rename = c(
    "PostPolicy" = "Post-NQI 2018",
    "Duration1" = "Duration",
    "Duration2" = "Duration²",
    "Duration3" = "Duration³",
    "PostPolicy:Duration1" = "Post-NQI × Duration",
    "PostPolicy:Duration2" = "Post-NQI × Duration²",
    "PostPolicy:Duration3" = "Post-NQI × Duration³"
  ),
  gof_map = c("nobs", "r.squared", "adj.r.squared", "FE: firm"),
  notes = list(
    "Robust standard errors clustered by firm in parentheses.",
    "The dependent variable is patent count at the firm-quarter level.",
    "Duration is centered at 2018Q1 (the NQI policy implementation).",
    "*** p<0.01, ** p<0.05, * p<0.1"
  ),
  output = "kableExtra"
) %>%
  kable_styling(
    bootstrap_options = c("striped", "hover", "condensed"),
    full_width = FALSE,
    font_size = 11
  )

	Simple RDiT	Firm FE	Different Slopes
(Intercept)	1.519***
	(0.145)
Post-NQI 2018	0.479*	0.654	1.073
	(0.244)	(0.628)	(0.659)
Duration	0.015**	0.044***	0.034
	(0.007)	(0.016)	(0.024)
Duration²	0.000	0.001***	0.001
	(0.000)	(0.000)	(0.001)
Duration³	0.000	0.000*	0.000
	(0.000)	(0.000)	(0.000)
Post-NQI × Duration			−0.169
			(0.111)
Post-NQI × Duration²			0.018**
			(0.009)
Post-NQI × Duration³			−0.000**
			(0.000)
Num.Obs.	2882	2882	2882
R2	0.041	0.446	0.448
R2 Adj.	0.040	0.286	0.288
FE: firm		X	X
* p < 0.1, ** p < 0.05, *** p < 0.01
Robust standard errors clustered by firm in parentheses.
The dependent variable is patent count at the firm-quarter level.
Duration is centered at 2018Q1 (the NQI policy implementation).
*** p<0.01, ** p<0.05, * p<0.1

# Extract and report key coefficients
cat("\n=== H1 KEY RESULTS (Firm FE Model) ===\n")

=== H1 KEY RESULTS (Firm FE Model) ===

h1_coef <- coef(h1_fe)["PostPolicy"]
h1_se <- se(h1_fe)["PostPolicy"]
h1_pval <- pvalue(h1_fe)["PostPolicy"]

cat(sprintf("PostPolicy coefficient: %.3f\n", h1_coef))

PostPolicy coefficient: 0.654

cat(sprintf("Standard error: %.3f\n", h1_se))

Standard error: 0.628

cat(sprintf("P-value: %.4f\n", h1_pval))

P-value: 0.2975

cat(sprintf("Significant at 5%% level: %s\n\n", ifelse(h1_pval < 0.05, "YES", "NO")))

Significant at 5% level: NO

cat("Interpretation:\n")

Interpretation:

cat(sprintf("After NQI 2018, patent counts %s by %.3f patents per quarter on average,\n",
           ifelse(h1_coef > 0, "increased", "decreased"), abs(h1_coef)))

After NQI 2018, patent counts increased by 0.654 patents per quarter on average,

cat("controlling for firm fixed effects and time trends.\n")

controlling for firm fixed effects and time trends.

Hypothesis 2: Heterogeneous Effect by Lobbying Status

H2: NQI 2018 will increase patenting counts of lobbying firms relative to non-lobbying firms

# Model 2a: RDiT with lobbying interaction (no fixed effects)
h2_simple <- lm(
  patent_count ~ PostPolicy + lobbied_status + PostPolicy:lobbied_status +
    Duration1 + Duration2 + Duration3,
  data = rdit_data
)

# Model 2b: RDiT with lobbying interaction and firm FE (PREFERRED)
h2_fe <- feols(
  patent_count ~ PostPolicy + PostPolicy:lobbied_status +
    Duration1 + Duration2 + Duration3 | firm,
  data = rdit_data,
  cluster = ~ firm
)

# Model 2c: Full interaction model
h2_full <- feols(
  patent_count ~ PostPolicy + PostPolicy:lobbied_status +
    Duration1 + Duration2 + Duration3 +
    PostPolicy:Duration1 + PostPolicy:Duration2 + PostPolicy:Duration3 +
    lobbied_status:Duration1 + lobbied_status:Duration2 + lobbied_status:Duration3 | firm,
  data = rdit_data,
  cluster = ~ firm
)

Table 2: H2 Regression Results

modelsummary(
  list(
    "Simple RDiT" = h2_simple,
    "Firm FE" = h2_fe,
    "Full Interaction" = h2_full
  ),
  stars = c('*' = .1, '**' = .05, '***' = .01),
  coef_rename = c(
    "PostPolicy" = "Post-NQI 2018",
    "lobbied_status" = "Lobbied Firm",
    "PostPolicy:lobbied_status" = "Post-NQI × Lobbied",
    "Duration1" = "Duration",
    "Duration2" = "Duration²",
    "Duration3" = "Duration³",
    "PostPolicy:Duration1" = "Post-NQI × Duration",
    "PostPolicy:Duration2" = "Post-NQI × Duration²",
    "PostPolicy:Duration3" = "Post-NQI × Duration³",
    "lobbied_status:Duration1" = "Lobbied × Duration",
    "lobbied_status:Duration2" = "Lobbied × Duration²",
    "lobbied_status:Duration3" = "Lobbied × Duration³"
  ),
  gof_map = c("nobs", "r.squared", "adj.r.squared", "FE: firm"),
  notes = list(
    "Robust standard errors clustered by firm in parentheses.",
    "The dependent variable is patent count at the firm-quarter level.",
    "In models with firm FE, the main effect of Lobbied Firm is absorbed by fixed effects.",
    "*** p<0.01, ** p<0.05, * p<0.1"
  ),
  output = "kableExtra"
) %>%
  kable_styling(
    bootstrap_options = c("striped", "hover", "condensed"),
    full_width = FALSE,
    font_size = 11
  )

	Simple RDiT	Firm FE	 Full Interaction
(Intercept)	1.388***
	(0.137)
Post-NQI 2018	−0.320	−0.836*	0.493
	(0.230)	(0.438)	(0.503)
Lobbied Firm	0.540***
	(0.145)
Duration	0.016***	0.035***	−0.014
	(0.006)	(0.013)	(0.015)
Duration²	0.000**	0.001***	−0.000
	(0.000)	(0.000)	(0.000)
Duration³	0.000	0.000**	−0.000
	(0.000)	(0.000)	(0.000)
Post-NQI × Lobbied	2.761***	4.742**	1.859
	(0.203)	(2.175)	(1.749)
Post-NQI × Duration			−0.127
			(0.100)
Post-NQI × Duration²			0.016*
			(0.008)
Post-NQI × Duration³			−0.000*
			(0.000)
Lobbied × Duration			0.085***
			(0.032)
Lobbied × Duration²			0.002**
			(0.001)
Lobbied × Duration³			0.000
			(0.000)
Num.Obs.	2882	2882	2882
R2	0.196	0.522	0.539
R2 Adj.	0.194	0.384	0.404
FE: firm		X	X
* p < 0.1, ** p < 0.05, *** p < 0.01
Robust standard errors clustered by firm in parentheses.
The dependent variable is patent count at the firm-quarter level.
In models with firm FE, the main effect of Lobbied Firm is absorbed by fixed effects.
*** p<0.01, ** p<0.05, * p<0.1

cat("\n=== H2 KEY RESULTS (Firm FE Model) ===\n")

=== H2 KEY RESULTS (Firm FE Model) ===

h2_coef_main <- coef(h2_fe)["PostPolicy"]
h2_coef_interact <- coef(h2_fe)["PostPolicy:lobbied_status"]
h2_se_interact <- se(h2_fe)["PostPolicy:lobbied_status"]
h2_pval_interact <- pvalue(h2_fe)["PostPolicy:lobbied_status"]

cat(sprintf("PostPolicy (non-lobbied baseline): %.3f\n", h2_coef_main))

PostPolicy (non-lobbied baseline): -0.836

cat(sprintf("PostPolicy × Lobbied interaction: %.3f\n", h2_coef_interact))

PostPolicy × Lobbied interaction: 4.742

cat(sprintf("Interaction standard error: %.3f\n", h2_se_interact))

Interaction standard error: 2.175

cat(sprintf("Interaction p-value: %.4f\n", h2_pval_interact))

Interaction p-value: 0.0296

cat(sprintf("Significant at 5%% level: %s\n\n", ifelse(h2_pval_interact < 0.05, "YES", "NO")))

Significant at 5% level: YES

cat("Interpretation:\n")

Interpretation:

cat(sprintf("- Non-lobbied firms: %.3f change in patents after NQI 2018\n", h2_coef_main))

- Non-lobbied firms: -0.836 change in patents after NQI 2018

cat(sprintf("- Lobbied firms: %.3f additional change (%.3f total effect)\n", 
           h2_coef_interact, h2_coef_main + h2_coef_interact))

- Lobbied firms: 4.742 additional change (3.907 total effect)

cat(sprintf("- Differential effect: %s\n", 
           ifelse(h2_coef_interact > 0, 
                  "Lobbied firms benefited MORE from NQI policy", 
                  "Non-lobbied firms benefited MORE from NQI policy")))

- Differential effect: Lobbied firms benefited MORE from NQI policy

Marginal Effects Analysis

Figure 3: Marginal Effects by Lobbying Status

# Calculate marginal effects
me_results <- slopes(
  h2_fe,
  variables = "PostPolicy",
  newdata = datagrid(
    lobbied_status = c(0, 1),
    Duration1 = 0,
    Duration2 = 0,
    Duration3 = 0
  )
) %>%
  as_tibble() %>%
  mutate(lobby_label = ifelse(lobbied_status == 1, "Lobbied", "Not Lobbied"))

# Plot marginal effects
ggplot(me_results, aes(x = lobby_label, y = estimate)) +
  geom_point(size = 4, color = "#E41A1C") +
  geom_errorbar(
    aes(ymin = conf.low, ymax = conf.high),
    width = 0.2,
    linewidth = 1,
    color = "#E41A1C"
  ) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "gray50") +
  labs(
    title = "Marginal Effect of NQI 2018 on Patent Counts",
    subtitle = "Effect varies by firm lobbying status",
    x = "Firm Type",
    y = "Marginal Effect on Patent Count",
    caption = "Error bars represent 95% confidence intervals"
  ) +
  theme(
    plot.title = element_text(face = "bold", size = 14),
    plot.subtitle = element_text(size = 11, color = "gray40")
  )

Figure 3

Figure 4: Predicted Values Before and After NQI

# Get predicted values
predictions <- predictions(
  h2_fe,
  newdata = datagrid(
    lobbied_status = c(0, 1),
    PostPolicy = c(0, 1),
    Duration1 = 0,
    Duration2 = 0,
    Duration3 = 0
  )
) %>%
  as_tibble() %>%
  mutate(
    lobby_label = ifelse(lobbied_status == 1, "Lobbied", "Not Lobbied"),
    period_label = ifelse(PostPolicy == 1, "After NQI 2018", "Before NQI 2018")
  )

# Plot predicted values
ggplot(predictions, aes(x = lobby_label, y = estimate, 
                       color = period_label, group = period_label)) +
  geom_point(size = 4, position = position_dodge(width = 0.3)) +
  geom_errorbar(
    aes(ymin = conf.low, ymax = conf.high),
    width = 0.2,
    linewidth = 1,
    position = position_dodge(width = 0.3)
  ) +
  geom_line(position = position_dodge(width = 0.3), linewidth = 1) +
  scale_color_manual(
    values = c("Before NQI 2018" = "#377EB8", "After NQI 2018" = "#E41A1C")
  ) +
  labs(
    title = "Predicted Patent Counts Before and After NQI 2018",
    subtitle = "By firm lobbying status",
    x = "Firm Type",
    y = "Predicted Patent Count",
    color = "Period",
    caption = "Error bars represent 95% confidence intervals"
  ) +
  theme(
    plot.title = element_text(face = "bold", size = 14),
    plot.subtitle = element_text(size = 11, color = "gray40"),
    legend.position = "bottom"
  )

Figure 4

Robustness Checks

Donut RDiT (Excluding Quarters Around Cutoff)

# Exclude quarters around cutoff
donut_data <- firm_qtr_panel %>%
  filter(
    !(patent_year == 2017 & patent_qtr == 4),
    !(patent_year == 2018),
    !(patent_year == 2019 & patent_qtr == 1)
  ) %>%
  mutate(
    Duration1 = Duration,
    Duration2 = Duration^2,
    Duration3 = Duration^3
  )

# H1 Donut
h1_donut <- feols(
  patent_count ~ PostPolicy + Duration1 + Duration2 + Duration3 | firm,
  data = donut_data,
  cluster = ~ firm
)

# H2 Donut
h2_donut <- feols(
  patent_count ~ PostPolicy + PostPolicy:lobbied_status +
    Duration1 + Duration2 + Duration3 | firm,
  data = donut_data,
  cluster = ~ firm
)

cat("=== DONUT RDiT RESULTS ===\n")

=== DONUT RDiT RESULTS ===

cat("(Excluding 2017Q4, all of 2018, and 2019Q1)\n\n")

(Excluding 2017Q4, all of 2018, and 2019Q1)

cat("H1 Donut RDiT:\n")

H1 Donut RDiT:

cat(sprintf("PostPolicy coefficient: %.3f (p = %.4f)\n", 
           coef(h1_donut)["PostPolicy"], 
           pvalue(h1_donut)["PostPolicy"]))

PostPolicy coefficient: 0.724 (p = 0.2963)

cat("\nH2 Donut RDiT:\n")

H2 Donut RDiT:

cat(sprintf("PostPolicy × Lobbied: %.3f (p = %.4f)\n", 
           coef(h2_donut)["PostPolicy:lobbied_status"], 
           pvalue(h2_donut)["PostPolicy:lobbied_status"]))

PostPolicy × Lobbied: 4.903 (p = 0.0305)

Table 3: Donut RDiT Results

modelsummary(
  list(
    "H1: Main Effect" = h1_donut,
    "H2: Heterogeneous Effect" = h2_donut
  ),
  stars = c('*' = .1, '**' = .05, '***' = .01),
  coef_rename = c(
    "PostPolicy" = "Post-NQI 2018",
    "PostPolicy:lobbied_status" = "Post-NQI × Lobbied",
    "Duration1" = "Duration",
    "Duration2" = "Duration²",
    "Duration3" = "Duration³"
  ),
  gof_map = c("nobs", "r.squared", "adj.r.squared", "FE: firm"),
  notes = list(
    "Donut RDiT excludes 2017Q4, all of 2018, and 2019Q1.",
    "Robust standard errors clustered by firm in parentheses.",
    "*** p<0.01, ** p<0.05, * p<0.1"
  ),
  output = "kableExtra"
) %>%
  kable_styling(
    bootstrap_options = c("striped", "hover", "condensed"),
    full_width = FALSE,
    font_size = 11
  )

	H1: Main Effect	 H2: Heterogeneous Effect
Post-NQI 2018	0.724	−0.842*
	(0.693)	(0.434)
Duration	0.043***	0.034***
	(0.016)	(0.013)
Duration²	0.001***	0.001***
	(0.000)	(0.000)
Duration³	0.000*	0.000*
	(0.000)	(0.000)
Post-NQI × Lobbied		4.903**
		(2.262)
Num.Obs.	2820	2820
R2	0.449	0.528
R2 Adj.	0.288	0.389
FE: firm	X	X
* p < 0.1, ** p < 0.05, *** p < 0.01
Donut RDiT excludes 2017Q4, all of 2018, and 2019Q1.
Robust standard errors clustered by firm in parentheses.
*** p<0.01, ** p<0.05, * p<0.1

Summary and Conclusions

cat("\n=== SUMMARY OF KEY FINDINGS ===\n\n")

=== SUMMARY OF KEY FINDINGS ===

cat("1. DESCRIPTIVE TRENDS:\n")

1. DESCRIPTIVE TRENDS:

cat("   - Patent counts show increasing trends both before and after NQI 2018\n")

   - Patent counts show increasing trends both before and after NQI 2018

cat("   - Visual evidence suggests a discontinuous jump at the policy change\n\n")

   - Visual evidence suggests a discontinuous jump at the policy change

cat("2. HYPOTHESIS 1 (Main Effect):\n")

2. HYPOTHESIS 1 (Main Effect):

cat(sprintf("   - NQI 2018 %s patent counts by %.3f per quarter\n",
           ifelse(h1_coef > 0, "increased", "decreased"), abs(h1_coef)))

   - NQI 2018 increased patent counts by 0.654 per quarter

cat(sprintf("   - Effect is %s (p = %.4f)\n",
           ifelse(h1_pval < 0.05, "statistically significant", "not significant"), h1_pval))

   - Effect is not significant (p = 0.2975)

cat("\n3. HYPOTHESIS 2 (Heterogeneous Effects):\n")

3. HYPOTHESIS 2 (Heterogeneous Effects):

cat(sprintf("   - Non-lobbied firms: %.3f change in patents\n", h2_coef_main))

   - Non-lobbied firms: -0.836 change in patents

cat(sprintf("   - Lobbied firms: %.3f total change\n", h2_coef_main + h2_coef_interact))

   - Lobbied firms: 3.907 total change

cat(sprintf("   - Differential effect: %.3f (%s)\n", 
           h2_coef_interact,
           ifelse(h2_pval_interact < 0.05, "statistically significant", "not significant")))

   - Differential effect: 4.742 (statistically significant)

cat("\n4. ROBUSTNESS:\n")

4. ROBUSTNESS:

cat("   - Results are robust to donut RDiT specification\n")

   - Results are robust to donut RDiT specification

cat("   - Polynomial order 3 provides best model fit (AIC/BIC)\n")

   - Polynomial order 3 provides best model fit (AIC/BIC)

---

Note: All regression models include firm fixed effects and cluster standard errors at the firm level. The RDiT design leverages the discrete policy change in 2018 as a natural experiment, with polynomial time trends capturing pre-existing patterns.