# Load necessary libraries
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.5
## ✔ forcats   1.0.0     ✔ stringr   1.5.1
## ✔ ggplot2   3.5.1     ✔ tibble    3.2.1
## ✔ lubridate 1.9.3     ✔ tidyr     1.3.1
## ✔ purrr     1.0.2     
## ── 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 errors
# Simulate data
set.seed(123) # for reproducibility
n <- 100 # Number of observations

# Generate random data for the variables
data <- tibble(
  prtystrA = rnorm(n, mean = 0.5, sd = 0.1),      # Party strength variable
  expendA = rnorm(n, mean = 300, sd = 50),        # Expenditure by candidate A
  expendB = rnorm(n, mean = 200, sd = 50),        # Expenditure by candidate B
  voteA = 50 + 10 * prtystrA + 0.05 * expendA - 0.03 * expendB + 
          0.001 * expendA * expendB + rnorm(n, sd = 5) # Vote share for candidate A
)

# View a sample of the data
head(data)
## # A tibble: 6 × 4
##   prtystrA expendA expendB voteA
##      <dbl>   <dbl>   <dbl> <dbl>
## 1    0.444    264.    310.  137.
## 2    0.477    313.    266.  142.
## 3    0.656    288.    187.  114.
## 4    0.507    283.    227.  121.
## 5    0.513    252.    179.  105.
## 6    0.672    298.    176.  120.
# Model with interaction between expendA and expendB
model <- lm(voteA ~ prtystrA + expendA + expendB + I(expendA * expendB), data = data)
summary(model)
## 
## Call:
## lm(formula = voteA ~ prtystrA + expendA + expendB + I(expendA * 
##     expendB), data = data)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -12.2272  -3.2248   0.2256   3.2429  12.3687 
## 
## Coefficients:
##                        Estimate Std. Error t value Pr(>|t|)    
## (Intercept)          73.6052762 14.8107962   4.970 2.96e-06 ***
## prtystrA              5.9540389  5.8472891   1.018   0.3111    
## expendA              -0.0203428  0.0473941  -0.429   0.6687    
## expendB              -0.1433059  0.0671279  -2.135   0.0353 *  
## I(expendA * expendB)  0.0013676  0.0002263   6.045 2.92e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 5.213 on 95 degrees of freedom
## Multiple R-squared:  0.9249, Adjusted R-squared:  0.9217 
## F-statistic: 292.4 on 4 and 95 DF,  p-value: < 2.2e-16
# Calculate average of expendA
avg_expendA <- mean(data$expendA, na.rm = TRUE)
print(paste("Average of expendA:", avg_expendA))
## [1] "Average of expendA: 294.622660099825"
# Fix expendA at 300 and estimate effect of additional $100,000 by Candidate B
expendA_fixed <- 300
increase_expendB <- 100
beta_3 <- coef(model)["expendB"]
beta_4 <- coef(model)["I(expendA * expendB)"]
predicted_change <- (beta_3 + beta_4 * expendA_fixed) * increase_expendB
print(paste("Effect of additional $100,000 by Candidate B on voteA:", predicted_change))
## [1] "Effect of additional $100,000 by Candidate B on voteA: 26.6988151537411"
# Fix expendB at 100 and examine effect of expendA = 100 on voteA
predicted_voteA <- predict(model, newdata = data.frame(prtystrA = 0, expendA = 100, expendB = 100))
print(paste("Predicted voteA when expendA is 100 and expendB is 100:", predicted_voteA))
## [1] "Predicted voteA when expendA is 100 and expendB is 100: 70.9168741904244"
# Replace interaction with shareA (Candidate A's share of total expenditures)
data <- data %>% mutate(shareA = expendA / (expendA + expendB))
model_shareA <- lm(voteA ~ prtystrA + expendA + shareA, data = data)
summary(model_shareA)
## 
## Call:
## lm(formula = voteA ~ prtystrA + expendA + shareA, data = data)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -14.6447  -3.8770   0.4805   3.3454  17.2698 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  120.72787    6.72133  17.962   <2e-16 ***
## prtystrA      11.73307    7.14915   1.641    0.104    
## expendA        0.42675    0.01599  26.688   <2e-16 ***
## shareA      -217.24925   11.75650 -18.479   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 6.424 on 96 degrees of freedom
## Multiple R-squared:  0.8847, Adjusted R-squared:  0.8811 
## F-statistic: 245.6 on 3 and 96 DF,  p-value: < 2.2e-16
# Partial effect of expendB on voteA with expendA = 300 and expendB = 0
partial_effect_expendB <- coef(model_shareA)["expendA"] - coef(model_shareA)["shareA"] * expendA_fixed / (expendA_fixed + 0)^2
print(paste("Partial effect of expendB on voteA with expendA = 300 and expendB = 0:", partial_effect_expendB))
## [1] "Partial effect of expendB on voteA with expendA = 300 and expendB = 0: 1.15091639541365"