rm(list=ls())
# Load necessary library
library(plm)
## Warning: 程辑包'plm'是用R版本4.2.3 来建造的
# Set seed for reproducibility
set.seed(123)

# Simulate a panel data set
N <- 100  # Number of individuals (cross-sectional units)
T <- 5    # Number of time periods

# Create a data frame with individual (ID) and time (Year) variables
data <- data.frame(
  ID = rep(1:N, each = T),
  Year = rep(1:T, times = N)
)

# Generate independent variables (X1, X2) and an error term (epsilon)
data$X1 <- rnorm(N * T, mean = 10, sd = 2)
data$X2 <- rnorm(N * T, mean = 5, sd = 1)
data$epsilon <- rnorm(N * T, mean = 0, sd = 1)

# Generate a dependent variable (Y) using a linear relationship
data$Y <- 3 + 1.5 * data$X1 + 0.5 * data$X2 + data$epsilon

# Convert the data frame to a panel data frame
pdata <- pdata.frame(data, index = c("ID", "Year"))

# Define a model formula
model_formula <- Y ~ X1 + X2

# Estimate a fixed effects model using the plm function
fixed_effects_model <- plm(model_formula, data = pdata, model = "within")

# Summary of the fixed effects model
summary(fixed_effects_model)
## Oneway (individual) effect Within Model
## 
## Call:
## plm(formula = model_formula, data = pdata, model = "within")
## 
## Balanced Panel: n = 100, T = 5, N = 500
## 
## Residuals:
##       Min.    1st Qu.     Median    3rd Qu.       Max. 
## -2.6384161 -0.6033977 -0.0055078  0.5788630  2.9062338 
## 
## Coefficients:
##    Estimate Std. Error t-value  Pr(>|t|)    
## X1 1.521652   0.024762  61.450 < 2.2e-16 ***
## X2 0.512405   0.048714  10.519 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Total Sum of Squares:    4126.3
## Residual Sum of Squares: 386.04
## R-Squared:      0.90645
## Adj. R-Squared: 0.8827
## F-statistic: 1928.1 on 2 and 398 DF, p-value: < 2.22e-16